← Back to branch summary

~ubuntu-branches/ubuntu/quantal/psicode/quantal

« back to all changes in this revision

Viewing changes to src/bin/ccresponse/transL.c

  • Committer: Bazaar Package Importer
  • Author(s): Michael Banck
  • Date: 2006-09-10 14:01:33 UTC
  • Revision ID: james.westby@ubuntu.com-20060910140133-ib2j86trekykfsfv
Tags: upstream-3.2.3
Import upstream version 3.2.3

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/bin/ccresponse/transL.c
 
1
#include <stdio.h>
 
2
#include <stdlib.h>
 
3
#include <libciomr/libciomr.h>
 
4
#include <libqt/qt.h>
 
5
#include <libiwl/iwl.h>
 
6
#define EXTERN
 
7
#include <psifiles.h>
 
8
#include "globals.h"
 
9
 
 
10
void transL(double sign)
 
11
{
 
12
  int stat, nao, noei_ao, nso, nmo;
 
13
  int i, j, ij, I, J;
 
14
  double *scratch, **TMP, **X;
 
15
  double **LX, **LY, **LZ;  /* MO-basis dipole integrals */
 
16
 
 
17
  nao = moinfo.nao;
 
18
  nso = moinfo.nso;
 
19
  nmo = moinfo.nmo;
 
20
  noei_ao = moinfo.noei_ao;
 
21
 
 
22
  /**** Transform the magnetic dipole integrals to the MO basis ****/
 
23
 
 
24
  LX = block_matrix(nmo,nmo);
 
25
  LY = block_matrix(nmo,nmo);
 
26
  LZ = block_matrix(nmo,nmo);
 
27
 
 
28
  TMP = block_matrix(nao, nao);
 
29
  X = block_matrix(nao, nao);
 
30
  scratch = init_array(noei_ao);
 
31
 
 
32
  /* NB: (a|m|b) = -1/2 (a|L|b) */
 
33
  /* NB: The angular momentum integrals are antisymmetric! */
 
34
  stat = iwl_rdone(PSIF_OEI, PSIF_AO_LX, scratch, noei_ao, 0, 0, outfile); /* read lower triangle */
 
35
  for(i=0,ij=0; i < nao; i++)
 
36
    for(j=0; j <= i; j++, ij++) {
 
37
      TMP[i][j] =  0.5 * scratch[ij] * sign;
 
38
      TMP[j][i] = -0.5 * scratch[ij] * sign;
 
39
    }
 
40
 
 
41
  C_DGEMM('n','t',nao,nso,nao,1,&(TMP[0][0]),nao,&(moinfo.usotao[0][0]),nao,
 
42
          0,&(X[0][0]),nao);
 
43
  C_DGEMM('n','n',nso,nso,nao,1,&(moinfo.usotao[0][0]),nao,&(X[0][0]),nao,
 
44
          0,&(TMP[0][0]),nao);
 
45
 
 
46
  C_DGEMM('n','n',nso,nmo,nso,1,&(TMP[0][0]),nao,&(moinfo.scf[0][0]),nmo,
 
47
          0,&(X[0][0]),nao);
 
48
  C_DGEMM('t','n',nmo,nmo,nso,1,&(moinfo.scf[0][0]),nmo,&(X[0][0]),nao,
 
49
          0,&(LX[0][0]),nmo);
 
50
 
 
51
  zero_arr(scratch,noei_ao);
 
52
 
 
53
  stat = iwl_rdone(PSIF_OEI, PSIF_AO_LY, scratch, noei_ao, 0, 0, outfile);
 
54
  for(i=0,ij=0; i < nao; i++)
 
55
    for(j=0; j <= i; j++, ij++) {
 
56
      TMP[i][j] =  0.5 * scratch[ij] * sign;
 
57
      TMP[j][i] = -0.5 * scratch[ij] * sign;
 
58
    }
 
59
 
 
60
  C_DGEMM('n','t',nao,nso,nao,1,&(TMP[0][0]),nao,&(moinfo.usotao[0][0]),nao,
 
61
          0,&(X[0][0]),nao);
 
62
  C_DGEMM('n','n',nso,nso,nao,1,&(moinfo.usotao[0][0]),nao,&(X[0][0]),nao,
 
63
          0,&(TMP[0][0]),nao);
 
64
 
 
65
  C_DGEMM('n','n',nso,nmo,nso,1,&(TMP[0][0]),nao,&(moinfo.scf[0][0]),nmo,
 
66
          0,&(X[0][0]),nao);
 
67
  C_DGEMM('t','n',nmo,nmo,nso,1,&(moinfo.scf[0][0]),nmo,&(X[0][0]),nao,
 
68
          0,&(LY[0][0]),nmo);
 
69
 
 
70
  zero_arr(scratch,noei_ao);
 
71
 
 
72
  stat = iwl_rdone(PSIF_OEI, PSIF_AO_LZ, scratch, noei_ao, 0, 0, outfile);
 
73
  for(i=0,ij=0; i < nao; i++)
 
74
    for(j=0; j <= i; j++, ij++) {
 
75
      TMP[i][j] =  0.5 * scratch[ij] * sign;
 
76
      TMP[j][i] = -0.5 * scratch[ij] * sign;
 
77
    }
 
78
 
 
79
  C_DGEMM('n','t',nao,nso,nao,1,&(TMP[0][0]),nao,&(moinfo.usotao[0][0]),nao,
 
80
          0,&(X[0][0]),nao);
 
81
  C_DGEMM('n','n',nso,nso,nao,1,&(moinfo.usotao[0][0]),nao,&(X[0][0]),nao,
 
82
          0,&(TMP[0][0]),nao);
 
83
 
 
84
  C_DGEMM('n','n',nso,nmo,nso,1,&(TMP[0][0]),nao,&(moinfo.scf[0][0]),nmo,
 
85
          0,&(X[0][0]),nao);
 
86
  C_DGEMM('t','n',nmo,nmo,nso,1,&(moinfo.scf[0][0]),nmo,&(X[0][0]),nao,
 
87
          0,&(LZ[0][0]),nmo);
 
88
 
 
89
  free_block(TMP);
 
90
  free_block(X);
 
91
  free(scratch);
 
92
 
 
93
  /*
 
94
  fprintf(outfile, "MO-Basis LX Integrals:\n");
 
95
  mat_print(LX, nmo, nmo, outfile);
 
96
  fprintf(outfile, "MO-Basis LY Integrals:\n");
 
97
  mat_print(LY, nmo, nmo, outfile);
 
98
  fprintf(outfile, "MO-Basis LZ Integrals:\n");
 
99
  mat_print(LZ, nmo, nmo, outfile);
 
100
 
 
101
  TMP = block_matrix(nmo, nmo);
 
102
  for(i=0; i < nmo; i++) {
 
103
    I = moinfo.pitzer2qt[i];
 
104
    for(j=0; j < nmo; j++) {
 
105
      J = moinfo.pitzer2qt[j];
 
106
      TMP[I][J] = LX[i][j];
 
107
    }
 
108
  }
 
109
  fprintf(outfile, "MO-Basis LX Integrals (QT):\n");
 
110
  print_mat(TMP, nmo, nmo, outfile);
 
111
  for(i=0; i < nmo; i++) {
 
112
    I = moinfo.pitzer2qt[i];
 
113
    for(j=0; j < nmo; j++) {
 
114
      J = moinfo.pitzer2qt[j];
 
115
      TMP[I][J] = LY[i][j];
 
116
    }
 
117
  }
 
118
  fprintf(outfile, "MO-Basis LY Integrals (QT):\n");
 
119
  print_mat(TMP, nmo, nmo, outfile);
 
120
  for(i=0; i < nmo; i++) {
 
121
    I = moinfo.pitzer2qt[i];
 
122
    for(j=0; j < nmo; j++) {
 
123
      J = moinfo.pitzer2qt[j];
 
124
      TMP[I][J] = LZ[i][j];
 
125
    }
 
126
  }
 
127
  fprintf(outfile, "MO-Basis LZ Integrals (QT):\n");
 
128
  print_mat(TMP, nmo, nmo, outfile);
 
129
  free_block(TMP);
 
130
  */
 
131
 
 
132
  moinfo.LX = LX;
 
133
  moinfo.LY = LY;
 
134
  moinfo.LZ = LZ;
 
135
}