1* GENCC production codes (whatever that means)
2      SUBROUTINE CI_TO_CC_REFRM(LUCC,LUCI,ISPC,ISM)
3*
4* A CI vector is given as the only vector in LUCI
5* Rewrite this vector to a set of Coupled CLuster amplitudes so
6*
7* Exp T |Ref> = CI
8*
9* Jeppe Olsen, April 14, early in the morning
10*
11* Reference space (CI space 1 ) is assumed to be a single det
12*
13* Output CC coefficients are put on FILE LUCC in current form
14* of CC coefficients.
15c      INCLUDE 'implicit.inc'
16c      INCLUDE 'mxpdim.inc'
17      INCLUDE 'wrkspc.inc'
18      INCLUDE 'clunit.inc'
19      INCLUDE 'cstate.inc'
20      INCLUDE 'csm.inc'
21      INCLUDE 'cicisp.inc'
22      INCLUDE 'glbbas.inc'
23      INCLUDE 'cgas.inc'
24      INCLUDE 'gasstr.inc'
25      INCLUDE 'strinp.inc'
26      INCLUDE 'crun.inc'
27      INCLUDE 'orbinp.inc'
28      INCLUDE 'ctcc.inc'
29      INCLUDE 'cprnt.inc'
30
31      CHARACTER*6 CCTYPE
32*
33      NTEST = 00
34*
35      IDUM = 0
36      CCTYPE(1:6) = 'GEN_CC'
37      CALL MEMMAN(IDUM,IDUM,'MARK  ',IDUM,'CICCRF')
38      CALL LUCIAQENTER('CI_CC ')
39*
40      WRITE(6,*) ' CI => CC transformation of coefficients '
41*. Space for CI behind the curtain
42      CALL GET_3BLKS_GCC(KVEC1,KVEC2,KVEC3,MXCJ)
43*
44* 1 : Find coefficient of reference det
45*
46      IREFSPC = 1
47      LBLK = -1
48C     EXPCIV(ISM,ISPCIN,LUIN,ISPCUT,LUUT,LBLK,
49C    &                  LUSCR,NROOT,ICOPY,IDC,NTESTG)
50      CALL EXPCIV(ISM,ISPC,LUCI,IREFSPC,LUSC1,LBLK,
51     &            LUHC,1,0,IDC,NTEST)
52      WRITE(6,*) ' LBLK after EXPCIV = ', LBLK
53*. Number of records in reference space
54      IF(IDC.EQ.1.OR.ISM.EQ.1) THEN
55        NREC = NSMST
56      ELSE
57        NREC = NSMST/2
58      END IF
59*. Read reference coefficient
60C          FRMDSCN(VEC,NREC,LBLK,LU)
61      CALL REWINO(LUSC1)
62      CALL FRMDSCN(CREF,NREC,LBLK,LUSC1)
63      WRITE(6,*) ' CREF = ', CREF
64      WRITE(6,*) ' LBLK after FRMDSCN = ', LBLK
65*
66* 2 : Normalize CI vector so reference coef is one
67*
68      IF(CREF.EQ.0.0D0) THEN
69        WRITE(6,*) ' CI_TO_CC_RF: Problems, norm of ref coef = 0'
70        STOP        'CI_TO_CC_RF: Problems, norm of ref coef = 0'
71      ELSE
72        FACTOR = 1.0D0/CREF
73*. And the scaling, result on LUHC
74        IREW = 1
75C            SCLVCD(LUIN,LUOUT,SCALE,SEGMNT,IREW,LBLK)
76        WRITE(6,*) ' Before call to SCLVCD '
77        WRITE(6,*) ' LUCI, LUHC = ', LUCI, LUHC
78        WRITE(6,*) ' KVEC1, LBLK = ', KVEC1, LBLK
79        CALL SCLVCD(LUCI,LUHC,FACTOR,WORK(KVEC1),IREW,LBLK)
80        WRITE(6,*) ' After call to SCLVCD '
81        IF(NTEST.GE.100) THEN
82          WRITE(6,*) ' Scaled CI vector '
83          CALL WRTVCD(WORK(KVEC1),LUHC,1,LBLK)
84        END IF
85      END IF
86*. So now : Complete vector with reference coef = 1 on LUHC
87*
88* 3. Information about CC space
89*
90*
91      IATP = 1
92      IBTP = 2
93*
94      NAEL = NELEC(IATP)
95      NBEL = NELEC(IBTP)
96      NEL = NAEL + NBEL
97*
98      IREFSPC = 1
99*. Find the type of reference state
100      CALL CC_AC_SPACES(IREFSPC,IREFTYP)
101*. Number of active orbital spaces
102      NACT_SPC = 0
103      IACT_SPC = 0
104      DO IGAS = 1, NGAS
105        IF(IHPVGAS(IGAS).EQ.3) THEN
106          NACT_SPC = NACT_SPC + 1
107          IACT_SPC = IGAS
108        END IF
109      END DO
110*. Info on active-active excitation types
111      CALL ACAC_EXC_TYP(IAAEXC_TYP,MX_AAEXC,IPRCC)
112*. Number of occupation classes for actual space
113      CALL OCCLSE(1,NOCCLS,IOCCLS,NEL,ISPC,0,0,NOBPT)
114*. And the occupation classes of actual space
115      CALL MEMMAN(KLOCCLS,NOCCLS*NEL,'ADDL  ',1,'OCCLS ')
116      CALL OCCLSE(2,NOCCLS,WORK(KLOCCLS),NEL,ISPC,0,0,NOBPT)
117*. Number of occupation classes for reference space
118      CALL OCCLSE(1,NOCCLS_REF,IOCCLS,NEL,IREFSPC,0,0,NOBPT)
119      IF(NOCCLS_REF.GT.1) THEN
120        WRITE(6,*) ' Problem in general CC '
121        WRITE(6,*)
122     &  ' Reference space is not a single occupation space'
123        STOP
124     &  ' Reference space is not a single occupation space'
125      END IF
126*. and the occupation classes of reference space
127      CALL MEMMAN(KLOCCLS_REF,NGAS,'ADDL  ',1,'OCC_RF')
128      CALL OCCLSE(2,NOCCLS_REF,WORK(KLOCCLS_REF),NEL,IREFSPC,0,0,NOBPT)
129*. Excitation type => Original occupation class
130*.
131*. Number of excitation types
132      IFLAG = 1
133      IDUM = 1
134      CALL TP_OBEX2(NOCCLS,NEL,NGAS,WORK(IDUM),
135     &             WORK(IDUM),WORK(IDUM),
136     &             WORK(KLOCCLS),WORK(KLOCCLS_REF),MX_NCREA,MX_NANNI,
137     &             MX_EXC_LEVEL,WORK(IDUM),MX_AAEXC,IFLAG,
138     &             I_OOCC,NOBEX_TP,NOAAEX,IPRCC)
139C?    WRITE(6,*) ' NOBEX_TP,MX_EXC_LEVEL = ', NOBEX_TP,MX_EXC_LEVEL
140*. And the actual orbital excitaions
141      NOBEX_TPE = NOBEX_TP + 1
142      CALL MEMMAN(KLCOBEX_TP,NOBEX_TPE,'ADDL  ',1,'LCOBEX')
143      CALL MEMMAN(KLAOBEX_TP,NOBEX_TPE,'ADDL  ',1,'LAOBEX')
144      CALL MEMMAN(KOBEX_TP ,NOBEX_TPE*2*NGAS,'ADDL  ',1,'IOBE_X')
145*. Excitation type => Original occupation class
146      CALL MEMMAN(KEX_TO_OC,NOBEX_TPE,'ADDL  ',1,'EX__OC')
147      IFLAG = 0
148      CALL TP_OBEX2(NOCCLS,NEL,NGAS,WORK(KOBEX_TP),
149     &             WORK(KLCOBEX_TP),WORK(KLAOBEX_TP),
150     &             WORK(KLOCCLS),WORK(KLOCCLS_REF),MX_NCREA,MX_NANNI,
151     &             MX_EXC_LEVEL,WORK(KEX_TO_OC),MX_AAEXC,IFLAG,
152     &             I_OOCC,NOBEX_TP,NOAAEX,IPRCC)
153*. Spinorbital excitations
154*. Spin combinations of CC excitations : Currently we assume that
155*. The T-operator is a singlet, can 'easily' be changed
156      IF(PSSIGN.EQ.0.0D0) THEN
157        MSCOMB_CC = 0
158      ELSE IF(PSSIGN.EQ.1.0D0) THEN
159        MSCOMB_CC = 1
160      END IF
161      MSCOMB_CC = 0
162*. Number of spin-orbital excitations
163      CALL OBEX_TO_SPOBEX(1,WORK(KOBEX_TP),WORK(KLCOBEX_TP),
164     &     WORK(KLAOBEX_TP),NOBEX_TP,IDUMMY,NSPOBEX_TP,NGAS,
165     &     NOBPT,0,MSCOMB_CC,IAAEXC_TYP,IACT_SPC,IPRCC,IDUMMY,
166     &     MXSPOX,WORK(KNSOX_FOR_OX),
167     &     WORK(KIBSOX_FOR_OX),WORK(KISOX_FOR_OX),
168     &     NAEL,NBEL,IREFSPC)
169*. And the actual spinorbital excitation operators
170      CALL MEMMAN(KLSOBEX,4*NGAS*NSPOBEX_TP,'ADDL  ',1,'SPOBEX')
171*. Map spin-orbital exc type => orbital exc type
172      CALL MEMMAN(KLSOX_TO_OX,NSPOBEX_TPE,'ADDL  ',1,'SPOBEX')
173*. First SOX of given OX ( including zero operator )
174      CALL MEMMAN(KIBSOX_FOR_OX,NOBEX_TP+1,'ADDL  ',1,'IBSOXF')
175*. Number of SOX's for given OX
176      CALL MEMMAN(KNSOX_FOR_OX,NOBEX_TP+1,'ADDL  ',1,'IBSOXF')
177*. SOX for given OX
178      CALL MEMMAN(KISOX_FOR_OX,NSPOBEX_TP+1,'ADDL  ',1,'IBSOXF')
179
180*. Map spin-orbital exc type => orbital exc type
181      CALL MEMMAN(KLSOX_TO_OX,NSPOBEX_TP,'ADDL  ',1,'SPOBEX')
182      CALL OBEX_TO_SPOBEX(2,WORK(KOBEX_TP),WORK(KLCOBEX_TP),
183     &     WORK(KLAOBEX_TP),NOBEX_TP,WORK(KLSOBEX),NSPOBEX_TP,NGAS,
184     &     NOBPT,0,MSCOMB_CC,IAAEXC_TYP,IACT_SPC,
185     &     IPRCC,WORK(KLSOX_TO_OX),MXSPOX,WORK(KNSOX_FOR_OX),
186     &     WORK(KIBSOX_FOR_OX),WORK(KISOX_FOR_OX),NAEL,NBEL,IREFSPC)
187*. Alpha- and beta-excitations constituting the spinorbital excitations
188*. Number
189      CALL SPOBEX_TO_ABOBEX(WORK(KLSOBEX),NSPOBEX_TP,NGAS,
190     &     1,NAOBEX_TP,NBOBEX_TP,IDUMMY,IDUMMY)
191*. And the alpha-and beta-excitations
192      LENA = 2*NGAS*NAOBEX_TP
193      LENB = 2*NGAS*NBOBEX_TP
194      CALL MEMMAN(KLAOBEX,LENA,'ADDL  ',2,'IAOBEX')
195      CALL MEMMAN(KLBOBEX,LENB,'ADDL  ',2,'IAOBEX')
196      CALL SPOBEX_TO_ABOBEX(WORK(KLSOBEX),NSPOBEX_TP,NGAS,
197     &     0,NAOBEX_TP,NBOBEX_TP,WORK(KLAOBEX),WORK(KLBOBEX))
198*. Max dimensions of CCOP !KSTR> = !ISTR> maps
199*. For alpha excitations
200      IATP = 1
201      IOCTPA = IBSPGPFTP(IATP)
202      NOCTPA = NSPGPFTP(IATP)
203      CALL LEN_GENOP_STR_MAP(
204     &     NAOBEX_TP,WORK(KLAOBEX),NOCTPA,NELFSPGP(1,IOCTPA),
205     &     NOBPT,NGAS,MAXLENA)
206      IBTP = 2
207      IOCTPB = IBSPGPFTP(IBTP)
208      NOCTPB = NSPGPFTP(IBTP)
209      CALL LEN_GENOP_STR_MAP(
210     &     NBOBEX_TP,WORK(KLBOBEX),NOCTPB,NELFSPGP(1,IOCTPB),
211     &     NOBPT,NGAS,MAXLENB)
212      MAXLEN_I1 = MAX(MAXLENA,MAXLENB)
213      WRITE(6,*) ' MAXLEN_I1 = ', MAXLEN_I1
214
215
216* Dimension of spinorbital excitation operators
217      ITOTSM = 1
218      CALL MEMMAN(KLLSOBEX,NSPOBEX_TP,'ADDL  ',1,'LSPOBX')
219      CALL MEMMAN(KLIBSOBEX,NSPOBEX_TP,'ADDL  ',1,'LSPOBX')
220      CALL MEMMAN(KLSPOBEX_AC,NSPOBEX_TP,'ADDL  ',1,'SPOBAC')
221*
222      CALL IDIM_TCC(WORK(KLSOBEX),NSPOBEX_TP,ITOTSM,
223     &          MX_ST_TSOSO_MX,MX_ST_TSOSO_BLK_MX,MX_TBLK_MX,
224     &          WORK(KLLSOBEX),WORK(KLIBSOBEX),LEN_T_VEC,
225     &          MSCOMB_CC,MX_SBSTR,
226     &          WORK(KISOX_FOR_OCCLS),NOCCLS,WORK(KIBSOX_FOR_OCCLS),
227     &          NTCONF,IPRCC)
228
229      N_CC_AMP = LEN_T_VEC
230      WRITE(6,*) 'N_CC_AMP = ', N_CC_AMP
231*. Allocate three CC vectors
232      CALL MEMMAN(KCCF,N_CC_AMP,'ADDL  ',2,'CCF   ')
233      CALL MEMMAN(KCC1,N_CC_AMP,'ADDL  ',2,'CC1   ')
234      CALL MEMMAN(KCC2,N_CC_AMP,'ADDL  ',2,'CC2   ')
235*
236      CALL MEMMAN(KLLCC,NSPOBEX_TP,'ADDL  ',1,'LCC   ')
237      CALL MEMMAN(KLICC,NSPOBEX_TP,'ADDL  ',1,'ICC   ')
238      CALL MEMMAN(KLJCC,NSPOBEX_TP,'ADDL  ',1,'JCC   ')
239*
240* Now the rest of the show goes as
241* Vector |LUHC> starts as complete CI with coef of refeence = 1
242*
243* Loop over excitation levels IEXC
244*  Reform CI vector !LUHC> to CC form
245*  Extract coefficients of excitation level IEXC,
246*  These are the CC coeffcients for this level
247*  Calculate Exp(-T(iexc))!LUHC> and store on LUHC
248* End of loop over excitation levels.
249*
250      ZERO = 0.0D0
251      CALL SETVEC(WORK(KCCF),ZERO,N_CC_AMP)
252      DO IEXC = 1, MX_EXC_LEVEL
253       IF(NTEST.GE.100) WRITE(6,*) ' Excitation level = ', IEXC
254*. Reform current CI coefficient to CC form, and store in WORK(KCC1)
255       IREW = 1
256       I_DO_CC_INFO = 0
257       CALL CC_CI_REORD(WORK(KCC1),LUHC,2,ISPC,ISM,IREW,I_DO_CC_INFO)
258*. Copy coefficients of excitation level to vector containing
259*. final CC amplitudes
260*. Spinorbital excitation types corresponding to this excitation level
261C     GET_SPOBTP_FOR_EXC_LEVEL(ILEVEL,ILEVEL_FOR_EXTP,
262C     &           NEXTP,NEXTP_AC,IEXTP_AC,ISOX_TO_OX)
263      CALL GET_SPOBTP_FOR_EXC_LEVEL(IEXC,WORK(KLCOBEX_TP),NSPOBEX_TP,
264     &     NEXTP_AC,WORK(KLJCC),WORK(KLSOX_TO_OX))
265*^ The active spinorbital excitation types are stored in WORK(KLJCC)
266*. first gathering from KCC1 to KCC2
267C     SCAGAT_CCVEC(CC_CMP,CC_EXP,ISG,NEXTP_SG,IEXTP_SG,
268C    &           IBEXTP,LEXTP,LEXTP_SG)
269       CALL  SCAGAT_CCVEC(WORK(KCC2),WORK(KCC1),2,NEXTP_AC,
270     &              WORK(KLJCC),WORK(KLIBSOBEX),WORK(KLLSOBEX),
271     &              WORK(KLLCC) )
272*. Then scatter from KCC1 to KCCF
273       CALL  SCAGAT_CCVEC(WORK(KCC2),WORK(KCCF),1,NEXTP_AC,
274     &              WORK(KLJCC),WORK(KLIBSOBEX),WORK(KLLSOBEX),
275     &              WORK(KLLCC) )
276*
277       IF(NTEST.GE.100) THEN
278         WRITE(6,*) ' Updated list of final CC coefficients : '
279         CALL WRTMAT(WORK(KCCF),1,N_CC_AMP,1,N_CC_AMP)
280       END IF
281*. calculate Exp(-T(iexc)|LUHC>
282*. Make only excitations with excitation level IEXC level
283       IZERO = 0
284       CALL ISETVC(WORK(KLSPOBEX_AC),IZERO,NSPOBEX_TP)
285       IONE = 1
286       CALL ISCASET(WORK(KLSPOBEX_AC),IONE,WORK(KLJCC),NEXTP_AC)
287       IF(NTEST.GE.100) THEN
288         WRITE(6,*) ' List of active spobex fresh from ISCASET '
289         CALL IWRTMA(WORK(KLSPOBEX_AC),1,NSPOBEX_TP,1,NSPOBEX_TP)
290       END IF
291*.  Exp(-t) !LUHC on LUSC35
292       MX_TERM = 100
293       ICC_EXC = 1
294       XCONV = 1.0D-20
295       CALL COPVEC(WORK(KCCF),WORK(KCC1),N_CC_AMP)
296       ONEM = -1.0D0
297       CALL SCALVE(WORK(KCC1),ONEM,N_CC_AMP)
298       CALL EXPT_REF(LUHC,LUSC35,LUSC1,LUSC2,LUSC3,XCONV,MX_TERM,
299     &               WORK(KVEC1),WORK(KVEC2),CCTYPE)
300*. And transfer to LUHC
301       CALL COPVCD(LUSC35,LUHC,WORK(KVEC1),1,LBLK)
302       IF(NTEST.GE.100) THEN
303         WRITE(6,*) ' Updated CI vector on LUHC '
304         CALL WRTVCD(WORK(KVEC1),LUHC,1,LBLK)
305       END IF
306      END DO
307*
308      IF(NTEST.GE.100) THEN
309       WRITE(6,*) ' CC coefficents obtained from CI coefficients'
310       CALL WRT_CC_VEC2(WORK(KCCF),IDUMMY,CCTYPE)
311      END IF
312*. Dump to LU_CI_TO_CC
313C     CALL REWINO(LU_CC_FROM_CI)
314C     WRITE(LU_CC_FROM_CI,'(I9)') N_CC_AMP
315C     DO I = 1, N_CC_AMP
316C       WRITE(LU_CC_FROM_CI,'(E25.15)') WORK(KCCF-1+I)
317C     END DO
318*. Dump to LUCCAMP
319      CALL REWINO(LU_CC_FROM_CI)
320      I_FORMATTED = 0
321      IF(I_FORMATTED.EQ.1) THEN
322        WRITE(LU_CC_FROM_CI,'(I9)') N_CC_AMP
323        DO I = 1, N_CC_AMP
324        WRITE(LU_CC_FROM_CI,'(E25.15)') WORK(KCCF-1+I)
325        END DO
326      ELSE
327        WRITE(LU_CC_FROM_CI) N_CC_AMP
328        WRITE(LU_CC_FROM_CI) (WORK(KCCF-1+I),I=1, N_CC_AMP)
329      END IF
330      CALL REWINO(LU_CC_FROM_CI)
331*
332      CALL MEMMAN(IDUM,IDUM,'FLUSM ',IDUM,'CICCRF')
333      CALL LUCIAQEXIT('CI_CC ')
334      RETURN
335      END
336      SUBROUTINE ISCASET(IARRAY,IVAL,ISCA,NSCA)
337*
338* IARRAY(ISCA(I)) = IVAL
339*
340* Jeppe Olsen, Aril 2000
341*
342      IMPLICIT REAL*8(A-H,O-Z)
343*. Input
344      INTEGER ISCA(NSCA)
345*. Output
346      INTEGER IARRAY(*)
347*
348      DO I = 1, NSCA
349       IARRAY(ISCA(I)) = IVAL
350      END DO
351*
352      RETURN
353      END
354      SUBROUTINE CC_CI_REORD(CCVEC,LUCI,IWAY,ISPC,ISM,IREW,
355     &           I_DO_CC_INFO)
356*
357* Convert between CI and CC organizations of coupled cluster
358* coefficients. Note that in this routine, the coefficients
359* are only reordered, no exponentations are involved here.
360*
361* Reference state is assumed to be a single Slaterdeterminant
362*
363* Input
364* =====
365*
366* CCVEC : Amplitudes organized as a CC vector
367* LUCI : File containing initial/final CI coefficients
368* IWAY  = 1 => CC to CI
369*       = 2 => CI to CC
370* ISPC : Space of expansions
371* ISM  : Symmetry of expansions
372*
373* Note : In core version in line with current assumption
374* that all coefs can be stored in core
375*
376* CI coefficients are initially/finally on file LUCI
377* but are in the routine store in a single array
378*
379* Jeppe Olsen, Magistratsvaegen 37, March 25  2000
380*              - in the kitchen, smelling Dittes cake and
381*                listening to Stones, Sticky  Fingers
382*
383c      INCLUDE 'implicit.inc'
384c      INCLUDE 'mxpdim.inc'
385      INCLUDE 'wrkspc.inc'
386      INCLUDE 'cgas.inc'
387      INCLUDE 'gasstr.inc'
388      INCLUDE 'orbinp.inc'
389      INCLUDE 'cstate.inc'
390      INCLUDE 'cicisp.inc'
391      INCLUDE 'strinp.inc'
392      INCLUDE 'stinf.inc'
393      INCLUDE 'strbas.inc'
394      INCLUDE 'csm.inc'
395      INCLUDE 'ctcc.inc'
396      INCLUDE 'crun.inc'
397*
398      DIMENSION CCVEC(*)
399*
400      IDUM = 0
401      CALL MEMMAN(IDUM,IDUM,'MARK  ',IDUM,'CC_CI_')
402*
403      NTEST = 00
404      LBLK = -1
405      IF(IREW.EQ.1) THEN
406        CALL REWINO(LUCI)
407      END IF
408*
409* Info on occupation classes in expansion
410*
411      IATP = 1
412      IBTP = 2
413*
414      NAEL = NELEC(IATP)
415      NBEL = NELEC(IBTP)
416      NEL = NAEL + NBEL
417*
418      ICSPC = ISPC
419      ISSPC = ISPC
420*
421      IREFSPC = 1
422*
423*. Number of occupation classes in CI and CC  expansion
424      CALL OCCLSE(1,NOCCLS,IOCCLS,NEL,ISPC,0,0,NOBPT)
425COLD  CALL OCCLS (1,NOCCLS,IOCCLS,NEL,NGAS,
426COLD &           IGSOCCX(1,1,ISPC),IGSOCCX(1,2,ISPC),
427COLD &           0,0,NOBPT)
428*. And the occupation classes
429      CALL MEMMAN(KLOCCLS,NOCCLS*NEL,'ADDL  ',1,'OCCLS ')
430      CALL OCCLSE(2,NOCCLS,WORK(KLOCCLS),NEL,ISPC,0,0,NOBPT)
431COLD  CALL OCCLS (2,NOCCLS,WORK(KLOCCLS),NEL,NGAS,
432COLD &           IGSOCCX(1,1,ISPC),IGSOCCX(1,2,ISPC),
433COLD &           0,0,NOBPT)
434*
435      IF(I_DO_CC_INFO.EQ.1) THEN
436*
437*              ==========================
438*              Information about CC space
439*              ==========================
440*
441*
442* Find the type of reference state
443*
444      CALL CC_AC_SPACES(IREFSPC,IREFTYP)
445*. Number of active orbital spaces
446      NACT_SPC = 0
447      IACT_SPC = 0
448      DO IGAS = 1, NGAS
449        IF(IHPVGAS(IGAS).EQ.3) THEN
450          NACT_SPC = NACT_SPC + 1
451          IACT_SPC = IGAS
452        END IF
453      END DO
454*. Info on active-active excitation types
455      CALL ACAC_EXC_TYP(IAAEXC_TYP,MX_AAEXC,IPRCC)
456*. Number of occupation classes for reference space
457      IREFSPC = 1
458      CALL OCCLSE(1,NOCCLS_REF,IOCCLS,NEL,IREFSPC,0,0,NOBPT)
459      IF(NOCCLS_REF.GT.1) THEN
460        WRITE(6,*) ' Problem in general CC '
461        WRITE(6,*)
462     &  ' Reference space is not a single occupation space'
463        STOP
464     &  ' Reference space is not a single occupation space'
465      END IF
466*. and the reference occupation space
467      CALL MEMMAN(KLOCCLS_REF,NGAS,'ADDL  ',1,'OCC_RF')
468      CALL OCCLSE(2,NOCCLS_REF,WORK(KLOCCLS_REF),NEL,IREFSPC,0,0,NOBPT)
469COLD  CALL OCCLS (2,NOCCLS_REF,WORK(KLOCCLS_REF),NEL,NGAS,
470COLD &           IGSOCCX(1,1,IREFSPC),IGSOCCX(1,2,IREFSPC),
471COLD &           0,0,NOBPT)
472*. Number of excitation types
473      IFLAG = 1
474      IDUM = 1
475      CALL TP_OBEX2(NOCCLS,NEL,NGAS,WORK(IDUM),
476     &             WORK(IDUM),WORK(IDUM),
477     &             WORK(KLOCCLS),WORK(KLOCCLS_REF),MX_NCREA,MX_NANNI,
478     &             MX_EXC_LEVEL,WORK(IDUM),MX_AAEXC,IFLAG,
479     &             I_OOCC,NOBEX_TP,NOAAEX,IPRCC)
480      WRITE(6,*) ' NOBEX_TP,MX_EXC_LEVEL = ', NOBEX_TP,MX_EXC_LEVEL
481      CALL MEMMAN(KLCOBEX_TP,NOBEX_TP,'ADDL  ',1,'LCOBEX')
482      CALL MEMMAN(KLAOBEX_TP,NOBEX_TP,'ADDL  ',1,'LAOBEX')
483      CALL MEMMAN(KOBEX_TP ,NOBEX_TP*2*NGAS,'ADDL  ',1,'IOBE_X')
484*. Excitation type => Original occupation class
485      CALL MEMMAN(KEX_TO_OC,NOBEX_TP,'ADDL  ',1,'EX__OC')
486      IFLAG = 0
487      CALL TP_OBEX2(NOCCLS,NEL,NGAS,WORK(KOBEX_TP),
488     &             WORK(KLCOBEX_TP),WORK(KLAOBEX_TP),
489     &             WORK(KLOCCLS),WORK(KLOCCLS_REF),MX_NCREA,MX_NANNI,
490     &             MX_EXC_LEVEL,WORK(KEX_TO_OC),MX_AAEXC,IFLAG,
491     &             I_OOCC,NOBEX_TP,NOAAEX,IPRCC)
492*. Spinorbital excitation types
493*. Spin combinations of CC excitations : Currently we assume that
494*. The T-operator is a singlet, can 'easily' be changed
495      IF(PSSIGN.EQ.0.0D0) THEN
496        MSCOMB_CC = 0
497      ELSE IF(PSSIGN.EQ.1.0D0) THEN
498        MSCOMB_CC = 1
499      END IF
500      MSCOMB_CC = 0
501*. Number of spinorbital excitation types
502      CALL OBEX_TO_SPOBEX(1,WORK(KOBEX_TP),WORK(KLCOBEX_TP),
503     &     WORK(KLAOBEX_TP),NOBEX_TP,IDUMMY,NSPOBEX_TP,NGAS,
504     &     NOBPT,0,MSCOMB_CC,IAAEXC_TYP,IACT_SPC,IPRCC,IDUMMY,
505     &     NAEL,NBEL)
506*. And the actual spinorbital excitation types
507      CALL MEMMAN(KLSOBEX,4*NGAS*NSPOBEX_TP,'ADDL  ',1,'SPOBEX')
508*. Map spin-orbital exc type => orbital exc type
509      CALL MEMMAN(KLSOX_TO_OX,NSPOBEX_TP,'ADDL  ',1,'SPOBEX')
510      CALL OBEX_TO_SPOBEX(2,WORK(KOBEX_TP),WORK(KLCOBEX_TP),
511     &     WORK(KLAOBEX_TP),NOBEX_TP,WORK(KLSOBEX),NSPOBEX_TP,NGAS,
512     &     NOBPT,0,MSCOMB_CC,IAAEXC_TYP,IACT_SPC,
513     &     IPRCC,WORK(KLSOX_TO_OX),NAEL,NBEL)
514* Dimension of spinorbital excitation operators
515      ITOTSM = 1
516      CALL MEMMAN(KLLSOBEX,NSPOBEX_TP,'ADDL  ',1,'LSPOBX')
517      CALL MEMMAN(KLIBSOBEX,NSPOBEX_TP,'ADDL  ',1,'LSPOBX')
518      CALL MEMMAN(KLSPOBEX_AC,NSPOBEX_TP,'ADDL  ',1,'SPOBAC')
519*
520      CALL IDIM_TCC(WORK(KLSOBEX),NSPOBEX_TP,ITOTSM,
521     &              MX_ST_TSOSO,MX_ST_TSOSO_BLK,MX_TBLK,
522     &              WORK(KLLSOBEX),WORK(KLIBSOBEX),LEN_T_VEC,
523     &              MSCOMB_CC,MX_SBSTR,IPRCC)
524      N_CC_AMP = LEN_T_VEC
525      END IF
526*
527*                     ==========================
528*                      Info for CI coefficients
529*                     ==========================
530*
531*
532*. Information about block structure- needed by new PICO2 routine.
533*. Memory for partitioning of C vector
534      NOCTPA = NOCTYP(IATP)
535      NOCTPB = NOCTYP(IBTP)
536      NTTS = MXNTTS
537C?    WRITE(6,*) ' GASCI : NTTS = ', NTTS
538      CALL MEMMAN(KLCLBT ,NTTS  ,'ADDL  ',1,'CLBT  ')
539      CALL MEMMAN(KLCLEBT ,NTTS  ,'ADDL  ',1,'CLEBT ')
540      CALL MEMMAN(KLCI1BT,NTTS  ,'ADDL  ',1,'CI1BT ')
541      CALL MEMMAN(KLCIBT ,8*NTTS,'ADDL  ',1,'CIBT  ')
542      CALL MEMMAN(KLC2B  ,  NTTS,'ADDL  ',1,'C2BT  ')
543*. Additional info required to construct partitioning
544*. Additional info required to construct partitioning
545*
546*
547      CALL MEMMAN(KLCIOIO,NOCTPA*NOCTPB,'ADDL  ',2,'CIOIO ')
548      CALL MEMMAN(KLCBLTP,NSMST,'ADDL  ',2,'CBLTP ')
549*
550      CALL IAIBCM(ISPC,WORK(KLCIOIO))
551*. option KSVST not active so
552      KSVST = 1
553      CALL ZBLTP(ISMOST(1,ISM),NSMST,IDC,WORK(KLCBLTP),WORK(KSVST))
554*
555*. Batches  of C vector
556      ICOMP = 1
557      ISIMSYM = 0
558*. Length of batch does not matter as we specified complete CI vector
559      LBLOCK = 1
560      CALL PART_CIV2(IDC,WORK(KLCBLTP),WORK(KNSTSO(IATP)),
561     &              WORK(KNSTSO(IBTP)),
562     &              NOCTPA,NOCTPB,NSMST,LBLOCK,WORK(KLCIOIO),
563     &              ISMOST(1,ISM),
564     &              NBATCH,WORK(KLCLBT),WORK(KLCLEBT),
565     &              WORK(KLCI1BT),WORK(KLCIBT),ICOMP,ISIMSYM)
566*. Number of BLOCKS
567        NBLOCK = IFRMR(WORK(KLCI1BT),1,NBATCH)
568     &         + IFRMR(WORK(KLCLBT),1,NBATCH) - 1
569C?      WRITE(6,*) ' Number of blocks ', NBLOCK
570*. Length of each block
571      CALL EXTRROW(WORK(KLCIBT),8,8,NBLOCK,WORK(KLCI1BT))
572*. Length of CI expansion
573      LENGTH_CI = IELSUM(WORK(KLCI1BT),NBLOCK)
574*. alphasupergroup, betasupergroup=> class
575        CALL MEMMAN(KLSPSPCL,NOCTPA*NOCTPB,'ADDL  ',1,'SPSPCL')
576        CALL SPSPCLS(WORK(KLSPSPCL),WORK(KLOCCLS),NOCCLS)
577*. Class of each block
578        CALL MEMMAN(KLBLKCLS,NBLOCK,'ADDL  ',1,'BLKCLS')
579        CALL MEMMAN(KLCLSL,NOCCLS,'ADDL  ',1,'CLSL  ')
580        CALL MEMMAN(KLCLSLR,NOCCLS,'ADDL  ',2,'CLSL_R  ')
581        CALL BLKCLS(WORK(KLCIBT),NBLOCK,WORK(KLBLKCLS),WORK(KLSPSPCL),
582     &              NOCCLS,WORK(KLCLSL),NOCTPA,NOCTPB,WORK(KLCLSLR))
583*
584* The connection between the CI and CC coefficients are
585* the mappings to the the Occupation classes
586*
587* KLBLKCLS : Occupation class for each CI block
588* KEXTP_TO_OCCLS : Occupation type for each excitation type
589*
590* Scratch vector for storing CI vector
591      CALL MEMMAN(KLCIVEC,LENGTH_CI,'ADDL  ',2,'CIVEC ')
592      IF(IWAY.EQ.2) THEN
593*. Collect CI vector in WORK(KLCIVEC)
594        IF(NTEST.GE.100) THEN
595          WRITE(6,*) ' Vector from LUCI '
596          CALL WRTVCD(WORK(KLCIVEC),LUCI,1,-1)
597          WRITE(6,*) ' LUCI,NBLOCK = ', LUCI,NBLOCK
598        END IF
599        CALL REWINO(LUCI)
600        LBLK = -1
601C       FRMDSCN(VEC,NREC,LBLK,LU)
602        CALL FRMDSCN(WORK(KLCIVEC),NBLOCK,LBLK,LUCI)
603      END IF
604*. Four blocks of string occupations for creation/annihilation strings
605      WRITE(6,*) ' MX_ST_TSOSO_BLK_MX = ', MX_ST_TSOSO_BLK_MX
606      CALL MEMMAN(KLSTR1_OCC,MX_ST_TSOSO_BLK_MX,'ADDL  ',1,'STOCC1')
607      CALL MEMMAN(KLSTR2_OCC,MX_ST_TSOSO_BLK_MX,'ADDL  ',1,'STOCC2')
608      CALL MEMMAN(KLSTR3_OCC,MX_ST_TSOSO_BLK_MX,'ADDL  ',1,'STOCC3')
609      CALL MEMMAN(KLSTR4_OCC,MX_ST_TSOSO_BLK_MX,'ADDL  ',1,'STOCC4')
610*. Space for string generation : Z matrices and strings
611*. Also used to hold an NORB*NORB matrix
612      LZSCR = (MAX(NAEL,NBEL)+3)*(NOCOB+1) + 2 * NOCOB + NOCOB*NOCOB
613      LZ    = (MAX(NAEL,NBEL)+2) * NOCOB
614      CALL MEMMAN(KLZSCR,LZSCR,'ADDL  ',2,'KLZSCR')
615      CALL MEMMAN(KLZ1,LZ,'ADDL  ',1,'KLZ1  ')
616      CALL MEMMAN(KLZ2,LZ,'ADDL  ',1,'KLZ2  ')
617*. Occupation af alpha- and betastrings
618      CALL MEMMAN(KLOCSTR1,MAX_STR_OC_BLK,'ADDL  ',1,'KLOCS1')
619      CALL MEMMAN(KLOCSTR2,MAX_STR_OC_BLK,'ADDL  ',1,'KLOCS2')
620*. Reorder arrays
621      CALL MEMMAN(KLREO1,MAX_STR_SPGP,'ADDL  ',1,'KLREO1')
622      CALL MEMMAN(KLREO2,MAX_STR_SPGP,'ADDL  ',1,'KLREO2')
623*. An alpha and betastring
624      CALL MEMMAN(KLSTRAL,NAEL,'ADDL  ',2,'STR_AL')
625      CALL MEMMAN(KLSTRBE,NBEL,'ADDL  ',2,'STR_BE')
626*
627      CALL CC_CI_REORD_S(CCVEC,WORK(KLCIVEC),IWAY,ISPC,ISM,
628     &     WORK(KLCIBT),NBLOCK,WORK(KLBLKCLS),
629     &     WORK(KLSOBEX),NSPOBEX_TP,WORK(KLSOX_TO_OX),
630     &     WORK(KLLSOBEX),WORK(KLIBSOBEX),WORK(KEX_TO_OC),
631     &     WORK(KLSTR1_OCC), WORK(KLSTR2_OCC), WORK(KLSTR3_OCC),
632     &     WORK(KLSTR4_OCC),WORK(KLZ1),WORK(KLZ2),
633     &     WORK(KLREO1),WORK(KLREO2),WORK(KLOCSTR1),WORK(KLOCSTR2),
634     &     WORK(KLZSCR),WORK(KLSTRAL),WORK(KLSTRBE),N_CC_AMP)
635
636      IF(IWAY.EQ.1) THEN
637*. Write resulting CI vector to DISC
638C       TODSCN(VEC,NREC,LREC,LBLK,LU)
639        CALL TODSCN(WORK(KLCIVEC),NBLOCK,WORK(KLCI1BT),LBLK,LUCI)
640      END IF
641*
642      IF(NTEST.GE.100) THEN
643        WRITE(6,*) ' Reordering between CI and CC order '
644        IF(IWAY.EQ.1) THEN
645          WRITE(6,*) ' CC to CI reordering '
646        ELSE IF ( IWAY.EQ.2) THEN
647          WRITE(6,*) ' CI to CC reordering '
648        END IF
649        WRITE(6,*) ' Vector of CC coefficients '
650        CALL WRTMAT(CCVEC,1,N_CC_AMP,1,N_CC_AMP)
651        WRITE(6,*) ' Vector of CI coefficients '
652        CALL WRTMAT(WORK(KLCIVEC),1,LENGTH_CI,1,LENGTH_CI)
653      END IF
654*
655COLD  STOP ' Enforced stop at end of CC_CI_REORD '
656      CALL MEMMAN(IDUM,IDUM,'FLUSM ',IDUM,'CC_CI_')
657      RETURN
658      END
659      SUBROUTINE SCAGAT_CCVEC(CC_CMP,CC_EXP,ISG,NEXTP_SG,IEXTP_SG,
660     &           IBEXTP,LEXTP,LEXTP_SG)
661*
662*. Scatter or gather blocks of CCVECTOR
663*  ISG = 1 : Scatter from ICC_CMP to ICC_EXP
664*  ISG = 2 : Gather  from ICC_EXP to ICC_CMP
665*
666* Blocks to be scattered/gathered are the NEXP_TP blocks in IEXP_TP
667*
668* Jeppe Olsen, April 2000
669*
670      INCLUDE 'implicit.inc'
671*. Input
672      INTEGER IEXTP_SG(NEXTP_SG)
673*. Input : Offset and length of spinorbital excitation blocks
674      INTEGER IBEXTP(*), LEXTP(*)
675*. Input/Output
676      DIMENSION CC_CMP(*),CC_EXP(*)
677*. Output (from gather)
678      INTEGER LEXTP_SG(*)
679*
680      NTEST = 00
681      IF(NTEST.GE.100) THEN
682        WRITE(6,*) 'SCAGAT_CCVEC speaking '
683      END IF
684
685      IOFF_CMP = 1
686      DO ITP = 1, NEXTP_SG
687        IOFF_EXP = IBEXTP(IEXTP_SG(ITP))
688        LEN = LEXTP(IEXTP_SG(ITP))
689        IF(NTEST.GE.100) THEN
690          WRITE(6,*) ' ITP, IOFF_EXP, LEN = ',ITP,IOFF_EXP,LEN
691        END IF
692        IF(ISG.EQ.1) THEN
693*. Scatter
694          CALL COPVEC(CC_CMP(IOFF_CMP),CC_EXP(IOFF_EXP),LEN)
695        ELSE
696*. Gather
697          CALL COPVEC(CC_EXP(IOFF_EXP),CC_CMP(IOFF_CMP),LEN)
698          LEXTP_SG(ITP) = LEN
699        END IF
700        IOFF_CMP = IOFF_CMP + LEN
701      END DO
702*
703      IF(NTEST.GE.100) THEN
704         WRITE(6,*) ' Gathered list, Vector and offsets '
705         LEN_CMP = IOFF_CMP-1
706         CALL WRTMAT(CC_CMP,1,LEN_CMP,1,LEN_CMP)
707         CALL IWRTMA(LEXTP_SG,1,NEXTP_SG,1,NEXTP_SG)
708      END IF
709*
710      RETURN
711      END
712      SUBROUTINE GET_SPOBTP_FOR_EXC_LEVEL(ILEVEL,ILEVEL_FOR_EXTP,
713     &           NEXTP,NEXTP_AC,IEXTP_AC,ISOX_TO_OX)
714*
715* Total number and blocknumbers of spinorbital excitations with
716* excitation level ILEVEL
717*
718* Jeppe Olsen, April 2000
719*
720      INCLUDE 'implicit.inc'
721*. Input
722      INTEGER ILEVEL_FOR_EXTP(NEXTP), ISOX_TO_OX(*)
723*. Output
724      INTEGER IEXTP_AC(*)
725*
726      NEXTP_AC = 0
727      DO IEXTP = 1, NEXTP
728*. Excitation level for this spinexcitation type
729        JLEVEL = ILEVEL_FOR_EXTP(ISOX_TO_OX(IEXTP))
730        IF(JLEVEL.EQ.ILEVEL) THEN
731          NEXTP_AC = NEXTP_AC + 1
732          IEXTP_AC(NEXTP_AC) = IEXTP
733        END IF
734      END DO
735*
736      NTEST = 00
737      IF(NTEST.GE.100) THEN
738        WRITE(6,*)
739     &  ' Spinorbital excitation blocks with excitation level',ILEVEL
740        WRITE(6,*)
741     &  ' Number of obtained spin-orbital excitation types', NEXTP_AC
742        WRITE(6,*) ' And the corresponding blocknumbers : '
743        CALL IWRTMA(IEXTP_AC,1,NEXTP_AC,1,NEXTP_AC)
744      END IF
745*
746      RETURN
747      END
748      SUBROUTINE CC_CI_REORD_S(CCVEC,CIVEC,IWAY,ISPC,ISM,
749     &           ICIBLK,NBLOCK_CI,ICIBLK_OCCLS,
750     &           ISOBEX,NSOBEX_TP,ISOX_TO_OX,LSOBEX,IBSOBEX,
751     &           IEX_TO_OC,
752     &           ICA_STR,ICB_STR,IAA_STR,IAB_STR,
753     &           IZA,IZB,IREOA,IREOB,IOCSTA,IOCSTB,IZSCR,
754     &           ISTRAL,ISTRBE,N_CC_AMP)
755*
756* Inner routine ( sounds nicer than slave routine )
757* for reordering between CI and CC orders.
758*
759* Only reordering is performed, no scaling
760*
761* Jeppe Olsen, March 28 2000
762*
763c      INCLUDE 'implicit.inc'
764c      INCLUDE 'mxpdim.inc'
765      INCLUDE 'wrkspc.inc'
766      INCLUDE 'cgas.inc'
767      INCLUDE 'gasstr.inc'
768      INCLUDE 'csm.inc'
769      INCLUDE 'multd2h.inc'
770      INCLUDE 'orbinp.inc'
771      INCLUDE 'strinp.inc'
772      INCLUDE 'newccp.inc'
773C     I_USE_NEWCCP
774*. Input and output
775      DIMENSION CCVEC(*),CIVEC(*)
776*. Input
777      INTEGER ICIBLK(8,NBLOCK_CI), ICIBLK_OCCLS(NBLOCK_CI)
778      INTEGER ISOBEX(4*NGAS,NSOBEX_TP),ISOX_TO_OX(NSOBEX_TP)
779      INTEGER LSOBEX(NSOBEX_TP),IBSOBEX(NSOBEX_TP)
780      INTEGER IEX_TO_OC(*)
781*. Space for creation and annihilation strings of given symmetry
782      INTEGER ICA_STR(*),ICB_STR(*),IAA_STR(*),IAB_STR(*)
783*. Space for strings, reorder arrays, and Z matrices, and scratch for
784*. constructing Z
785      INTEGER IZA(*),IZB(*),IREOA(*),IREOB(*),IOCSTA(*),IOCSTB(*)
786      INTEGER IZSCR(*)
787*. Space for a single alpha and beta string
788      INTEGER ISTRAL(*), ISTRBE(*)
789*. Local scratch : Occupation in Reference space
790      INTEGER IREF_OCC_AL(MXPNGAS),IREF_OCC_BE(MXPNGAS)
791*. Actual reference strings
792      INTEGER IREF_STR_AL(MXPNGAS),IREF_STR_BE(MXPNGAS)
793*. General occupation of a pair of alpha- and beta-strings
794      INTEGER IOCC_AL(MXPNGAS),IOCC_BE(MXPNGAS)
795*. And the corresponding groups
796      INTEGER IGRP_AL(MXPNGAS),IGRP_BE(MXPNGAS)
797*. Offsets to CI blocks, with given TT as a function of symmetry of
798*. alpha strings
799      INTEGER ICIBLK_OFF(MXPOBS)
800*. For group notation of annihilation/creation strings
801      INTEGER IGRP_CA(MXPNGAS),IGRP_CB(MXPNGAS)
802      INTEGER IGRP_AA(MXPNGAS),IGRP_AB(MXPNGAS)
803*. For local testing
804CTEST INTEGER ITOUCH(1000)
805*
806CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
807CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
808CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
809CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
810CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
811CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
812CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
813CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
814CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
815CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
816CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
817CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
818CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
819CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
820CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
821CTEST WRITE(6,*) ' Jeppe : Remember local tests are active '
822CTEST IZERO = 0
823CTEST CALL ISETVC(ITOUCH,IZERO,N_CC_AMP+1)
824*
825
826      ITP_AL = 1
827      ITP_BE = 2
828      NEL_AL = NELEC(ITP_AL)
829      NEL_BE = NELEC(ITP_BE)
830      I_CC = 0
831*
832      NTEST = 1000
833*. Check sums for CI and CC adressing
834      ICC_CHECK = 0
835      ICI_CHECK = 0
836*
837C?    WRITE(6,*) ' Included Spinorbital excitations'
838C?    CALL WRT_SPOX_TP(ISOBEX,NSOBEX_TP)
839      IDUM = 0
840      CALL MEMMAN(IDUM,IDUM,'MARK  ',IDUM,'CI_CIS')
841C     GET_REF_ALBE_OCC(IREFSPC,IREF_AL,IREF_BE)
842*. Obtain alpha and beta occupation of reference space
843      IREFSPC = 1
844      CALL GET_REF_ALBE_OCC(IREFSPC,IREF_OCC_AL,IREF_OCC_BE)
845*. Obtain atual alpha and beta strings for reference space
846      CALL GET_REF_ALBE_STR(IREFSPC,IREF_STR_AL,IREF_STR_BE)
847*. Symmetry of reference strings
848      ISM_REF_AL = ISYMST(IREF_STR_AL,NEL_AL)
849      ISM_REF_BE = ISYMST(IREF_STR_BE,NEL_BE)
850*. Loop over spinorbital excitation blocks
851      DO JSOBEX = 1, NSOBEX_TP
852C?     WRITE(6,*) ' Output for JSOBEX = ', JSOBEX
853*. Resulting occupation of alpha and beta strings
854C     EXOCC_STROCC(ICR_OCC,IAN_OCC,ISTR_IN_OCC,
855C    &           ISTR_OUT_OCC,NGAS,IZERO_STR)
856*. Occupation of alpha string
857        CALL EXOCC_STROCC(ISOBEX(1+0*NGAS,JSOBEX),
858     &       ISOBEX(1+2*NGAS,JSOBEX),IREF_OCC_AL,
859     &       IOCC_AL,NGAS,IZERO_ALSTR)
860        CALL OCC_TO_GRP(IOCC_AL,IGRP_AL,1)
861*. Occupation of betastring
862        CALL EXOCC_STROCC(ISOBEX(1+1*NGAS,JSOBEX),
863     &       ISOBEX(1+3*NGAS,JSOBEX),IREF_OCC_BE,
864     &       IOCC_BE,NGAS,IZERO_BESTR)
865        CALL OCC_TO_GRP(IOCC_BE,IGRP_BE,1)
866        IF(NTEST.GE.100) THEN
867          WRITE(6,*) ' Occupation of resulting strings for JSOBEX=',
868     &    JSOBEX
869          CALL IWRTMA(IOCC_AL,1,NGAS,1,NGAS)
870          CALL IWRTMA(IOCC_BE,1,NGAS,1,NGAS)
871          WRITE(6,*) ' And the corresponding groups '
872          CALL IWRTMA(IGRP_AL,1,NGAS,1,NGAS)
873          CALL IWRTMA(IGRP_BE,1,NGAS,1,NGAS)
874        END IF
875*. Supergroups corresponding to these occupation
876        CALL FIND_SPGRP_FROM_OCC(IOCC_AL,ISPGRP_AL,ITP_AL)
877        CALL FIND_SPGRP_FROM_OCC(IOCC_BE,ISPGRP_BE,ITP_BE)
878*. Relative number of these supergroups
879        ISPGRP_AL_REL = ISPGRP_AL - IBSPGPFTP(ITP_AL) + 1
880        ISPGRP_BE_REL = ISPGRP_BE - IBSPGPFTP(ITP_BE) + 1
881*. And then the TTS blocks with these occupation
882C            CIBLOCKS_FOR_TT(ICIBLK,NCIBLK,IATP,IBTP,IFORM,ITTBLK)
883        CALL CIBLOCKS_FOR_TT(ICIBLK,NBLOCK_CI,ISPGRP_AL_REL,
884     &       ISPGRP_BE_REL,2,ICIBLK_OFF)
885*. Transform creation/annihilations type from occupation to group notation
886        CALL OCC_TO_GRP(ISOBEX(1+0*NGAS,JSOBEX),IGRP_CA,1)
887        CALL OCC_TO_GRP(ISOBEX(1+1*NGAS,JSOBEX),IGRP_CB,1)
888        CALL OCC_TO_GRP(ISOBEX(1+2*NGAS,JSOBEX),IGRP_AA,1)
889        CALL OCC_TO_GRP(ISOBEX(1+3*NGAS,JSOBEX),IGRP_AB,1)
890*
891        NEL_CA = IELSUM(ISOBEX(1+0*NGAS,JSOBEX),NGAS)
892        NEL_CB = IELSUM(ISOBEX(1+1*NGAS,JSOBEX),NGAS)
893        NEL_AA = IELSUM(ISOBEX(1+2*NGAS,JSOBEX),NGAS)
894        NEL_AB = IELSUM(ISOBEX(1+3*NGAS,JSOBEX),NGAS)
895        IF(NTEST.GE.100) THEN
896          WRITE(6,*) ' NEL_CA, NEL_CB, NEL_AA, NEL_AB',
897     &                 NEL_CA, NEL_CB, NEL_AA, NEL_AB
898        END IF
899
900*. Loop over symmetries of creation/annihilation strings
901*. Symmetry of excitations is assumed to be 1 (total sym)
902        ISM = 1
903        DO ISM_C = 1, NSMST
904          ISM_A = MULTD2H(ISM,ISM_C)
905          DO ISM_CA = 1, NSMST
906            ISM_CB = MULTD2H(ISM_C,ISM_CA)
907            DO ISM_AA = 1, NSMST
908             ISM_AB =  MULTD2H(ISM_A,ISM_AA)
909*. Obtain creation/annihilation strings
910             CALL GETSTR2_TOTSM_SPGP(IGRP_CA,NGAS,ISM_CA,NEL_CA,
911     &            NSTR_CA,ICA_STR, NTOOB,0,IDUM,IDUM)
912             CALL GETSTR2_TOTSM_SPGP(IGRP_CB,NGAS,ISM_CB,NEL_CB,
913     &            NSTR_CB,ICB_STR, NTOOB,0,IDUM,IDUM)
914             CALL GETSTR2_TOTSM_SPGP(IGRP_AA,NGAS,ISM_AA,NEL_AA,
915     &            NSTR_AA,IAA_STR, NTOOB,0,IDUM,IDUM)
916             CALL GETSTR2_TOTSM_SPGP(IGRP_AB,NGAS,ISM_AB,NEL_AB,
917     &            NSTR_AB,IAB_STR, NTOOB,0,IDUM,IDUM)
918*
919C?       WRITE(6,*) ' Beta annihilation strings : '
920C?       CALL IWRTMA(IAB_STR,NEL_AB,NSTR_AB,NEL_AB,NSTR_AB)
921*. Corresponding symmetries of alpha and beta strings
922             ISM_OP_AL = MULTD2H(ISM_CA,ISM_AA)
923             ISM_OP_BE = MULTD2H(ISM_CB,ISM_AB)
924*. Symmetry of alpha and beta strings
925             ISM_STR_AL = MULTD2H(ISM_REF_AL,ISM_OP_AL)
926             ISM_STR_BE = MULTD2H(ISM_REF_BE,ISM_OP_BE)
927*. Obtain all alpha and beta strings of correct supergroup and sym
928*. The mapping from occupation to address will be used on the following
929*.. Generate information about IA strings
930C                 WEIGHT_SPGP(Z,NORBTP,NELFTP,NORBFTP,ISCR,NTEST)
931             NTESTL = 0
932             CALL WEIGHT_SPGP(IZA,NGAS,IOCC_AL,NOBPT,IZSCR,NTESTL)
933C            GETSTR2_TOTSM_SPGP(IGRP,NIGRP,ISPGRPSM,NEL,NSTR,ISTR,
934C    &                              NORBT,IDOREO,IZ,IREO)
935             CALL GETSTR2_TOTSM_SPGP(IGRP_AL,NGAS,ISM_STR_AL,NEL_AL,
936     &            NSTR_AL,IOCSTA,NOCOB,1,IZA,IREOA)
937C?           WRITE(6,*) ' Reorder array for alpha strings '
938C?           CALL IWRTMA(IREOA,1,NSTR_AL,1,NSTR_AL)
939*. And about beta string
940             CALL WEIGHT_SPGP(IZB,NGAS,IOCC_BE,NOBPT,IZSCR,NTESTL)
941             CALL GETSTR2_TOTSM_SPGP(IGRP_BE,NGAS,ISM_STR_BE,NEL_BE,
942     &            NSTR_BE,IOCSTB,NOCOB,1,IZB,IREOB)
943*. Loop over T elements as  matric T(I_CA, I_CB, IAA, I_AB)
944             DO I_AB = 1, NSTR_AB
945              DO I_AA = 1, NSTR_AA
946               DO I_CB = 1, NSTR_CB
947                DO I_CA = 1, NSTR_CA
948*. Alpha string obtained by alpha crea alpha anni alpha refstring
949                 IOFF_CA = 1 + (I_CA-1)*NEL_CA
950                 IOFF_AA = 1 + (I_AA-1)*NEL_AA
951                 CALL CRAN_STR(ICA_STR(IOFF_CA),IAA_STR(IOFF_AA),
952     &                NEL_CA,NEL_AA,IREF_STR_AL,NEL_AL,
953     &                ISTRAL,ISIGN_AL,IZERO_AL)
954*. And number of this string
955C       ISTRNM(IOCC,NORB,NEL,Z,NEWORD,IREORD)
956                 IANUM = ISTRNM(ISTRAL,NOCOB,NEL_AL,IZA,IREOA,1)
957                 IF(NTEST.GE.1000) WRITE(6,*) ' IANUM = ', IANUM
958*. And for beta string
959                 IOFF_CB = 1 + (I_CB-1)*NEL_CB
960                 IOFF_AB = 1 + (I_AB-1)*NEL_AB
961                 IF(NTEST.GE.1000) THEN
962                   WRITE(6,*) ' I_AB, IOFF_AB = ',
963     &                          I_AB, IOFF_AB
964                 END IF
965                 CALL CRAN_STR(ICB_STR(IOFF_CB),IAB_STR(IOFF_AB),
966     &                NEL_CB,NEL_AB,IREF_STR_BE,NEL_BE,
967     &                ISTRBE,ISIGN_BE,IZERO_BE)
968*. And number of this string
969                 IBNUM = ISTRNM(ISTRBE,NOCOB,NEL_BE,IZB,IREOB,1)
970                 IF(NTEST.GE.1000) WRITE(6,*) ' IBNUM = ', IBNUM
971
972*. Number in CC order
973                 I_CC = I_CC + 1
974                 ICC_CHECK = ICC_CHECK + I_CC
975*. Number in CI order
976                 I_CI = ICIBLK_OFF(ISM_STR_AL)-1+(IBNUM-1)*NSTR_AL+
977     &                  IANUM
978CTEST            ITOUCH(I_CI) = ITOUCH(I_CI) + 1
979                 ICI_CHECK = ICI_CHECK + I_CI
980                 IF(NTEST.GE.1000) THEN
981                   WRITE(6,*) 'ICIBLK_OFF, NSTR_AL = ',
982     &                          ICIBLK_OFF(ISM_STR_AL),NSTR_AL
983                   WRITE(6,'(A,4I4)') ' I_AB, I_AA, I_CB, I_CA',
984     &                                  I_AB, I_AA, I_CB, I_CA
985                   WRITE(6,*) 'I_CC, I_CI = ',I_CC,I_CI
986
987                 END IF
988*
989                 IF(I_USE_NEWCCP.EQ.0) THEN
990                   SIGN = DFLOAT(ISIGN_AL*ISIGN_BE)
991                 ELSE
992                   IF(MOD(NEL_CB*NEL_AA,2).EQ.0) THEN
993                     SIGN_CBAA = 1
994                   ELSE
995                     SIGN_CBAA = -1
996                   END IF
997                   SIGN = SIGN_CBAA* DFLOAT(ISIGN_AL*ISIGN_BE)
998                 END IF
999*
1000                 IF(IWAY.EQ.1) THEN
1001                   CIVEC(I_CI) = CCVEC(I_CC)*SIGN
1002                 ELSE
1003                   CCVEC(I_CC) = CIVEC(I_CI)*SIGN
1004                 END IF
1005*
1006                 IF(I_CI.LE.0.OR.I_CI.GT.N_CC_AMP+1) THEN
1007                   WRITE(6,*) ' I_CI out of range = ', I_CI
1008                 END IF
1009*
1010                END DO
1011               END DO
1012             END DO
1013*       C    ^ End of loop over elements of block
1014            END DO
1015           END DO
1016          END DO
1017        END DO
1018*       ^ End of loop over symmetries of creation/annihilation strings
1019      END DO
1020*     ^ End of loop over types of CC excitations
1021*
1022      IF(NTEST.GE.100) THEN
1023       WRITE(6,*) ' CC and CI check sums = ', ICC_CHECK,ICI_CHECK
1024      END IF
1025*
1026*. Local tests : Print the numbers for the CI coefficients that
1027*. were not touched exactly once
1028CTEST WRITE(6,*) ' Local tests : '
1029CTEST WRITE(6,*) ' Local tests : '
1030CTEST WRITE(6,*) ' Local tests : '
1031CTEST WRITE(6,*) ' Local tests : '
1032CTEST WRITE(6,*) ' Local tests : '
1033CTEST WRITE(6,*) ' Local tests : '
1034CTEST WRITE(6,*) ' Local tests : '
1035CTEST WRITE(6,*) ' Local tests : '
1036CTEST WRITE(6,*) ' Local tests : '
1037CTEST WRITE(6,*) ' Local tests : '
1038CTEST DO ICI = 1, N_CC_AMP + 1
1039CTEST   IF(ITOUCH(ICI).NE.1) THEN
1040CTEST     WRITE(6,*) ' Mess : ICI, ITOUCH(ICI) = ',
1041CTEST&                        ICI, ITOUCH(ICI)
1042CTEST   END IF
1043CTEST END DO
1044*
1045*. The part below does not work as we do not know the number of
1046*. the CI det corresponding to the HF reference
1047C     IF(ICI_CHECK.NE.ICC_CHECK+1) THEN
1048C       WRITE(6,*) ' Problem in reord, inconsistent checksums'
1049C       WRITE(6,*) ' ICC_CHECK, ICI_CHECK = ', ICC_CHECK,ICI_CHECK
1050C       STOP       ' Problem in reord, inconsistent checksums'
1051C     END IF
1052*
1053      CALL MEMMAN(IDUM,IDUM,'FLUSM ',IDUM,'CI_CIS')
1054      RETURN
1055      END
1056      SUBROUTINE GET_REF_ALBE_OCC(IREFSPC,IREF_AL,IREF_BE)
1057*
1058* Obtain alpha and beta occupations for reference space
1059*
1060* Reference space is assumed to be a single pair of occupations of
1061* alpha and beta strings ( this includes closed shell HF,
1062* Highspin open shell and CAS reference)
1063*
1064* Only a single valence orbital space is assumed
1065*
1066* Jeppe Olsen, March 2000
1067*
1068      INCLUDE 'implicit.inc'
1069      INCLUDE 'mxpdim.inc'
1070      INCLUDE 'cgas.inc'
1071      INCLUDE 'gasstr.inc'
1072      INCLUDE 'strinp.inc'
1073*. Output : Alpha and beta occupations for each GAS space
1074      INTEGER IREF_AL(NGAS),IREF_BE(NGAS)
1075*
1076*. Total number of Hole orbitals
1077      NHOLE = 0
1078      DO IGAS = 1, NGAS
1079        IF(IHPVGAS(IGAS).EQ.1) THEN
1080          NHOLE = NHOLE + NGSOBT(IGAS)
1081        END IF
1082      END DO
1083*. Number of orbitals in valence space
1084      NELEC_AL = NELEC(1)
1085      NELEC_BE = NELEC(2)
1086      NVAL_AL = NELEC_AL - NHOLE
1087      NVAL_BE = NELEC_BE - NHOLE
1088*
1089      DO IGAS = 1, NGAS
1090        NORB = NGSOBT(IGAS)
1091        IF(IHPVGAS(IGAS).EQ.1) THEN
1092          IREF_AL(IGAS) = NORB
1093          IREF_BE(IGAS) = NORB
1094        ELSE IF( IHPVGAS(IGAS).EQ.2) THEN
1095          IREF_AL(IGAS) = 0
1096          IREF_BE(IGAS) = 0
1097        ELSE IF( IHPVGAS(IGAS).EQ.3) THEN
1098          IREF_AL(IGAS) = NVAL_AL
1099          IREF_BE(IGAS) = NVAL_BE
1100        END IF
1101      END DO
1102*
1103      NTEST = 00
1104      IF(NTEST.GE.100) THEN
1105        WRITE(6,*) ' Alpha and Beta Occupation for reference space '
1106        WRITE(6,*)
1107        CALL IWRTMA(IREF_AL,1,NGAS,1,NGAS)
1108        CALL IWRTMA(IREF_BE,1,NGAS,1,NGAS)
1109      END IF
1110*
1111      RETURN
1112      END
1113      SUBROUTINE GET_REF_ALBE_STR(IREFSPC,IREF_AL,IREF_BE)
1114*
1115* Obtain alpha and beta strings for reference space
1116*
1117* Reference space is assumed to be a single pair of STRINGS of
1118* alpha and beta strings ( this includes closed shell HF,
1119* Highspin open shell, but nor CAS !!)
1120*
1121* Only a single valence orbital space is assumed
1122*
1123* Jeppe Olsen, March 2000
1124*
1125      INCLUDE 'implicit.inc'
1126      INCLUDE 'mxpdim.inc'
1127      INCLUDE 'cgas.inc'
1128      INCLUDE 'gasstr.inc'
1129      INCLUDE 'strinp.inc'
1130*. Output : Alpha and beta occupations for each GAS space
1131      INTEGER IREF_AL(NGAS),IREF_BE(NGAS)
1132*
1133*. Total number of Hole orbitals
1134      NHOLE = 0
1135      DO IGAS = 1, NGAS
1136        IF(IHPVGAS(IGAS).EQ.1) THEN
1137          NHOLE = NHOLE + NGSOBT(IGAS)
1138        END IF
1139      END DO
1140*. Number of orbitals in valence space
1141      NELEC_AL = NELEC(1)
1142      NELEC_BE = NELEC(2)
1143      NVAL_AL = NELEC_AL - NHOLE
1144      NVAL_BE = NELEC_BE - NHOLE
1145*
1146      IOFF = 1
1147      IOFF_AL = 1
1148      IOFF_BE = 1
1149*
1150      DO IGAS = 1, NGAS
1151        IF(IGAS.EQ.1) THEN
1152          IOFF = 1
1153        ELSE
1154          IOFF = IOFF + NGSOBT(IGAS-1)
1155        END IF
1156        NORB = NGSOBT(IGAS)
1157        IF(IHPVGAS(IGAS).EQ.1) THEN
1158          DO IORB = 1, NORB
1159            IREF_AL(IOFF_AL-1+IORB) = IOFF-1+IORB
1160            IREF_BE(IOFF_BE-1+IORB) = IOFF-1+IORB
1161          END DO
1162          IOFF_AL = IOFF_AL + NORB
1163          IOFF_BE = IOFF_BE + NORB
1164        ELSE IF( IHPVGAS(IGAS).EQ.3) THEN
1165          DO IORB = 1, NVAL_AL
1166            IREF_AL(IOFF_AL-1+IORB) = IOFF-1+IORB
1167          END DO
1168          IOFF_AL = IOFF_AL + NVAL_AL
1169          DO IORB = 1, NVAL_BE
1170            IREF_BE(IOFF_BE-1+IORB) = IOFF-1+IORB
1171          END DO
1172          IOFF_BE = IOFF_BE + NVAL_BE
1173        END IF
1174      END DO
1175*
1176      NTEST = 000
1177      IF(NTEST.GE.100) THEN
1178        WRITE(6,*) ' Alpha and Beta strings for reference space '
1179        WRITE(6,*)
1180        CALL IWRTMA(IREF_AL,1,NELEC_AL,1,NELEC_AL)
1181        CALL IWRTMA(IREF_BE,1,NELEC_BE,1,NELEC_BE)
1182      END IF
1183*
1184      RETURN
1185      END
1186      SUBROUTINE EXOCC_STROCC(ICR_OCC,IAN_OCC,ISTR_IN_OCC,
1187     &           ISTR_OUT_OCC,NGAS,IZERO_STR)
1188*
1189* Occupation of excitaion op,  occupation of string =>
1190* Occupation of excited string
1191*
1192* Jeppe Olsen, March 2000
1193*
1194      INCLUDE 'implicit.inc'
1195*. Input
1196      INTEGER ICR_OCC(*),IAN_OCC(*),ISTR_IN_OCC(*)
1197*. Output
1198      INTEGER ISTR_OUT_OCC(*)
1199*. Annihilation
1200      IZERO_STR = 0
1201      DO IGAS = 1, NGAS
1202        ISTR_OUT_OCC(IGAS) = ISTR_IN_OCC(IGAS) - IAN_OCC(IGAS)
1203        IF(ISTR_OUT_OCC(IGAS).LT.0) IZERO_STR = 1
1204      END DO
1205*. Creation
1206      DO IGAS = 1, NGAS
1207        ISTR_OUT_OCC(IGAS) = ISTR_OUT_OCC(IGAS) + ICR_OCC(IGAS)
1208      END DO
1209*
1210      NTEST = 00
1211      IF(NTEST.GE.100) THEN
1212        WRITE(6,*) ' Output from EXOCC_STROCC '
1213        WRITE(6,*) ' ========================='
1214        WRITE(6,*)
1215        WRITE(6,*) ' Occ of crea string : '
1216        CALL IWRTMA(ICR_OCC,1,NGAS,1,NGAS)
1217        WRITE(6,*) ' Occ of anni string '
1218        CALL IWRTMA(IAN_OCC,1,NGAS,1,NGAS)
1219        WRITE(6,*) ' Occ of input string '
1220        CALL IWRTMA(ISTR_IN_OCC,1,NGAS,1,NGAS)
1221        WRITE(6,*) ' Occ of output string '
1222        CALL IWRTMA(ISTR_OUT_OCC,1,NGAS,1,NGAS)
1223      END IF
1224*
1225      RETURN
1226      END
1227      SUBROUTINE CRAN_STR(ICR,IAN,NCR,NAN,ISTR_IN,NEL_IN,
1228     &                    ISTR_OUT,ISIGN,IZERO_STR)
1229*
1230* ISTR_OUT = ISIGN* ICR IAN ISTR_IN
1231*
1232* Where ICR is a string of creation operators and IAN is a string
1233* of annihilation operators.
1234*
1235* Input string is assumed to be given in ascending order,
1236* and output string will be delivered with orbitals in
1237* ascending order
1238*
1239*. Initial version, I hope it is not for mission critical routines
1240* (could be speeded up)
1241*
1242* Jeppe Olsen, March 2000
1243*
1244* Change of phase of annihilations strings, Oct2000
1245*
1246      INCLUDE 'implicit.inc'
1247      INCLUDE 'newccp.inc'
1248*. Input
1249      INTEGER ICR(NCR),IAN(NAN)
1250      INTEGER ISTR_IN(NEL_IN)
1251*. Output
1252      INTEGER ISTR_OUT(*)
1253*
1254      NTEST = 00
1255      IF(NTEST.GE.100) THEN
1256        WRITE(6,*) ' CRAN_STR speaking '
1257        WRITE(6,*) ' =================='
1258        WRITE(6,*) ' Input string '
1259        CALL IWRTMA(ISTR_IN,1,NEL_IN,1,NEL_IN)
1260        WRITE(6,*) ' Annihilation string '
1261        CALL IWRTMA(IAN,1,NAN,1,NAN)
1262        WRITE(6,*) ' Creation string '
1263        CALL IWRTMA(ICR,1,NCR,1,NCR)
1264      END IF
1265*. Make sure that annihilation strings are properly increasing
1266C?    DO JAN = 2, NAN
1267C?      IF(IAN(JAN).LE.IAN(JAN-1)) THEN
1268C?        WRITE(6,*) ' CRAN confused, strange annihilation string :'
1269C?        CALL IWRTMA(IAN,1,NAN,1,NAN)
1270C?        STOP ' CRAN confused, strange annihilation string'
1271C?      END IF
1272C?    END DO
1273*
1274      NEL_OUT = NEL_IN - NAN + NCR
1275*
1276      IZERO_STR = 0
1277      ISIGN = 1.0D0
1278      CALL ICOPVE(ISTR_IN,ISTR_OUT,NEL_IN)
1279*. Annihilate  : IAN(1) IAN(2) .... !STR_IN>
1280      DO IANNI = 1, NAN
1281        IFOUND = 0
1282        DO IEL = 1, NEL_IN-IANNI+1
1283C?        WRITE(6,*) ' CRAN : IANNI IEL ISTR IAN ',
1284C?   &    IANNI,IEL,ISTR_OUT(IEL),IAN(NAN-IANNI+1)
1285          IF(ISTR_OUT(IEL).EQ.IAN(NAN-IANNI+1)) THEN
1286            ISIGN = ISIGN*(-1)**(IEL-1)
1287            IFOUND = 1
1288            DO JEL = IEL, NEL_IN-IANNI
1289              ISTR_OUT(JEL) = ISTR_OUT(JEL+1)
1290             END DO
1291          END IF
1292        END DO
1293        IF(IFOUND.EQ.0) THEN
1294*. orbital to be annihilated not found, output string is zero
1295          IZERO_STR = 1
1296          GOTO 1001
1297        END  IF
1298      END DO
1299      IF(NTEST.GE.100) THEN
1300        WRITE(6,*) ' Annihilated string '
1301        CALL IWRTMA(ISTR_OUT,NEL_IN-NAN,1,NEL_IN-NAN,1)
1302      END IF
1303*. Creation part
1304      DO ICREA = 1, NCR
1305*. Place to insert orbital
1306        ICR_ORB = ICR(NCR-ICREA+1)
1307        IPLACE = 1
1308        DO IEL = 1, NEL_IN-NAN + ICREA-1
1309          IF(ISTR_OUT(IEL).EQ.ICR_ORB) THEN
1310*. Electron is already around, zero
1311            IZERO_STR = 1
1312            GOTO 1001
1313          END IF
1314*
1315          IF(IEL.LT.NEL_IN-NAN + ICREA-1) THEN
1316            IF(ISTR_OUT(IEL).LT.ICR_ORB.AND.
1317     &         ISTR_OUT(IEL+1).GT.ICR_ORB   ) THEN
1318               IPLACE = IEL+1
1319            END IF
1320          ELSE IF(IEL.EQ.NEL_IN-NAN+ICREA-1) THEN
1321            IF(ISTR_OUT(IEL).LT. ICR_ORB   ) THEN
1322              IPLACE = IEL + 1
1323            END IF
1324          END IF
1325        END DO
1326        ISIGN  = ISIGN*(-1)**(IPLACE-1)
1327        DO IEL = NEL_IN-NAN+ICREA,IPLACE + 1, -1
1328          ISTR_OUT(IEL) = ISTR_OUT(IEL-1)
1329        END DO
1330        ISTR_OUT(IPLACE) = ICR_ORB
1331      END DO
1332*
1333 1001 CONTINUE
1334*. A bit on the sign : In LUCIA the order of the annihilation is actually
1335* descending !, find permutation sign
1336      IF(I_USE_NEWCCP.EQ.0.AND.NAN.GT.1) THEN
1337        NPERM = NAN*(NAN-1)/2
1338        ISIGN = ISIGN*(-1)**NPERM
1339      END IF
1340*
1341      IF(NTEST.GE.100) THEN
1342        IF(IZERO_STR.EQ.0) THEN
1343          WRITE(6,*) ' Output string '
1344          CALL IWRTMA(ISTR_OUT,1,NEL_OUT,1,NEL_OUT)
1345          WRITE(6,*) ' ISIGN = ', ISIGN
1346        ELSE
1347          WRITE(6,*) ' Vanishing string '
1348        END IF
1349      END IF
1350*
1351      RETURN
1352      END
1353C
1354      SUBROUTINE FIND_SPGRP_FROM_OCC(IOCC,ISPGRP_NUM,ITP)
1355*
1356* Find the number(ISPGRP_NUM) corresponding to supergroup of type ITP
1357* with occupation IOCC. If ITP = 0, all supergroup types are checked
1358* Returned supergroup number is absolute supergroup number
1359*
1360* If no supergroup is identified a zero is returned
1361*
1362* Jeppe Olsen, April 2000
1363* ITP = 0 option added, March 2007 - not tested..
1364*
1365      INCLUDE 'implicit.inc'
1366*. General input
1367      INCLUDE 'mxpdim.inc'
1368      INCLUDE 'cgas.inc'
1369      INCLUDE 'gasstr.inc'
1370*. Specific input
1371      INTEGER IOCC(*)
1372*
1373      ISPGRP_NUM = 0
1374      IF(ITP.EQ.0) THEN
1375        ITP_MIN = 1
1376        ITP_MAX = NTSPGP
1377      ELSE
1378        ITP_MIN = ITP
1379        ITP_MAX = ITP
1380      END IF
1381*
1382      DO IITP = ITP_MIN, ITP_MAX
1383       DO  JSPGP = IBSPGPFTP(IITP), IBSPGPFTP(IITP)+NSPGPFTP(IITP)-1
1384         IDENTICAL = 1
1385         DO IGAS = 1, NGAS
1386           IF(NELFSPGP(IGAS,JSPGP).NE.IOCC(IGAS)) IDENTICAL = 0
1387         END DO
1388         IF(IDENTICAL.EQ.1) ISPGRP_NUM = JSPGP
1389       END DO
1390      END DO
1391*
1392      NTEST = 00
1393      IF(NTEST.GE.100) THEN
1394        WRITE(6,*) ' Occupation of supergroup : '
1395        CALL IWRTMA(IOCC,1,NGAS,1,NGAS)
1396        IF(ISPGRP_NUM.EQ.0) THEN
1397          WRITE(6,*) ' Not identified '
1398        ELSE
1399          WRITE(6,*) ' Number of supergroup : ', ISPGRP_NUM
1400        END IF
1401      END IF
1402*
1403      RETURN
1404      END
1405      SUBROUTINE CIBLOCKS_FOR_TT(ICIBLK,NCIBLK,IATP,IBTP,IFORM,ITTBLK)
1406*
1407* A set of CI blocks is specified through ICIBLK
1408* Find block with TYPES IATP, IBTP
1409*
1410*
1411* Output :  ITTBLK : The blocks : IFORM = 1 => The number of the block
1412*                                       = 2 => Offset of block
1413*
1414*. Jeppe Olsen, April 20000
1415*
1416*. General input
1417      INCLUDE 'implicit.inc'
1418      INCLUDE 'mxpdim.inc'
1419      INCLUDE 'csm.inc'
1420      INTEGER ICIBLK(8,NCIBLK)
1421*. Output
1422      INTEGER ITTBLK(*)
1423*
1424      IZERO = 0
1425      CALL ISETVC(ITTBLK,IZERO,NSMST)
1426      DO JCIBLK = 1, NCIBLK
1427C?      WRITE(6,*) ' JCIBLK, IA and IB : ',
1428C?   &  JCIBLK,ICIBLK(1,JCIBLK),ICIBLK(2,JCIBLK)
1429        IF(ICIBLK(1,JCIBLK).EQ.IATP.AND.ICIBLK(2,JCIBLK).EQ.IBTP) THEN
1430          IASM = ICIBLK(3,JCIBLK)
1431          IOFF = ICIBLK(5,JCIBLK)
1432          IF(IFORM.EQ.1) THEN
1433            ITTBLK(IASM) = JCIBLK
1434          ELSE IF(IFORM.EQ.2) THEN
1435            ITTBLK(IASM) = IOFF
1436          END IF
1437        END IF
1438      END DO
1439*
1440      NTEST = 00
1441      IF(NTEST.GE.100) THEN
1442        WRITE(6,'(A,2I6)') ' Blocks with IATP, IBTP ', IATP,IBTP
1443        IF(IFORM.EQ.1) THEN
1444          WRITE(6,*) ' Block numbers '
1445        ELSE
1446          WRITE(6,*) ' Offsets '
1447        END IF
1448        CALL IWRTMA(ITTBLK,1,NSMST,1,NSMST)
1449      END IF
1450*
1451      RETURN
1452      END
1453
1454c $Id$
1455