← Back to branch summary

~ubuntu-branches/ubuntu/vivid/psicode/vivid

« back to all changes in this revision

Viewing changes to src/bin/ccdensity/transdip.c

  • Committer: Bazaar Package Importer
  • Author(s): Michael Banck, Michael Banck, Daniel Leidert
  • Date: 2009-02-23 00:12:02 UTC
  • mfrom: (1.1.2 upstream)
  • Revision ID: james.westby@ubuntu.com-20090223001202-rutldoy3dimfpesc
Tags: 3.4.0-1
* New upstream release.

[ Michael Banck ]
* debian/patches/01_DESTDIR.dpatch: Refreshed.
* debian/patches/02_FHS.dpatch: Removed, applied upstream.
* debian/patches/03_debian_docdir: Likewise.
* debian/patches/04_man.dpatch: Likewise.
* debian/patches/06_466828_fix_gcc_43_ftbfs.dpatch: Likewise.
* debian/patches/07_464867_move_executables: Fixed and refreshed.
* debian/patches/00list: Adjusted.
* debian/control: Improved description.
* debian/patches-held: Removed.
* debian/rules (install/psi3): Do not ship the ruby bindings for now.

[ Daniel Leidert ]
* debian/rules: Fix txtdir via DEB_MAKE_INSTALL_TARGET.
* debian/patches/01_DESTDIR.dpatch: Refreshed.

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/bin/ccdensity/transdip.c
1
 
#include <stdio.h>
2
 
#include <stdlib.h>
3
 
#include <libciomr/libciomr.h>
4
 
#include <libiwl/iwl.h>
5
 
#include <libchkpt/chkpt.h>
6
 
#include <libdpd/dpd.h>
7
 
#include <libqt/qt.h>
8
 
#include <psifiles.h>
9
 
#define EXTERN
10
 
#include "globals.h"
11
 
 
12
 
#define IOFF_MAX 32641
13
 
#define INDEX(i,j) ((i>j) ? (ioff[(i)]+(j)) : (ioff[(j)]+(i)))
14
 
 
15
 
void transdip(void)
16
 
{
17
 
  int nmo, nso, nao, noei, stat, i, I, h, j;
18
 
  int *order, *doccpi, *ioff;
19
 
  double **scf_pitzer, **scfA, **scfB;
20
 
  double **scf_qt, **X, **usotao;
21
 
  double *zvals, **geom;
22
 
  double *mu_x_ints, *mu_y_ints, *mu_z_ints;
23
 
  double **MUX_AO, **MUY_AO, **MUZ_AO;
24
 
  double **MUX_MO, **MUY_MO, **MUZ_MO;
25
 
 
26
 
  chkpt_init(PSIO_OPEN_OLD);
27
 
  if ((params.ref == 0) || (params.ref == 1))
28
 
    scf_pitzer = chkpt_rd_scf();
29
 
  else if(params.ref == 2) {
30
 
    scfA = chkpt_rd_alpha_scf();
31
 
    scfB = chkpt_rd_beta_scf();
32
 
  }
33
 
 
34
 
  nao = chkpt_rd_nao();
35
 
  nso = chkpt_rd_nso();
36
 
  usotao = chkpt_rd_usotao();
37
 
  chkpt_close();
38
 
 
39
 
  nmo = moinfo.nmo;
40
 
 
41
 
  /*** Build ioff ***/
42
 
  ioff = init_int_array(IOFF_MAX);
43
 
  ioff[0] = 0;
44
 
  for(i=1; i < IOFF_MAX; i++) ioff[i] = ioff[i-1] + i;
45
 
 
46
 
  /*** Get the Pitzer -> QT reordering array ***/
47
 
  order = init_int_array(nmo);
48
 
 
49
 
  /* doccpi array must include frozen orbitals for reorder_qt() */
50
 
  doccpi = init_int_array(moinfo.nirreps);
51
 
  for(h=0; h < moinfo.nirreps; h++) 
52
 
    doccpi[h] = moinfo.frdocc[h] + moinfo.clsdpi[h];
53
 
 
54
 
  reorder_qt(doccpi, moinfo.openpi, moinfo.frdocc, moinfo.fruocc, 
55
 
             order, moinfo.orbspi, moinfo.nirreps);
56
 
 
57
 
  /*** Reorder the SCF eigenvectors to QT ordering */
58
 
 
59
 
  scf_qt = block_matrix(nmo, nmo);
60
 
  for(i=0; i < nmo; i++) {
61
 
    I = order[i];  /* Pitzer --> QT */
62
 
    for(j=0; j < nmo; j++) scf_qt[j][I] = scf_pitzer[j][i];
63
 
  }
64
 
  free(order);
65
 
  free(doccpi);
66
 
  free_block(scf_pitzer);
67
 
 
68
 
  /*** Read in dipole moment integrals in the AO basis ***/
69
 
  noei = nao * (nao + 1)/2;
70
 
 
71
 
  mu_x_ints = init_array(noei);
72
 
  stat = iwl_rdone(PSIF_OEI,PSIF_AO_MX,mu_x_ints,noei,0,0,outfile);
73
 
  mu_y_ints = init_array(noei);
74
 
  stat = iwl_rdone(PSIF_OEI,PSIF_AO_MY,mu_y_ints,noei,0,0,outfile);
75
 
  mu_z_ints = init_array(noei);
76
 
  stat = iwl_rdone(PSIF_OEI,PSIF_AO_MZ,mu_z_ints,noei,0,0,outfile);
77
 
 
78
 
  MUX_AO = block_matrix(nao,nao);
79
 
  MUY_AO = block_matrix(nao,nao);
80
 
  MUZ_AO = block_matrix(nao,nao);
81
 
 
82
 
  for(i=0; i < nao; i++)
83
 
    for(j=0; j < nao; j++) {
84
 
      MUX_AO[i][j] = mu_x_ints[INDEX(i,j)];
85
 
      MUY_AO[i][j] = mu_y_ints[INDEX(i,j)];
86
 
      MUZ_AO[i][j] = mu_z_ints[INDEX(i,j)];
87
 
    }
88
 
 
89
 
/*   fprintf(outfile, "MUX_AOs\n"); */
90
 
/*   print_mat(MUX_AO, nao, nao, outfile); */
91
 
/*   fprintf(outfile, "MUY_AOs\n"); */
92
 
/*   print_mat(MUY_AO, nao, nao, outfile); */
93
 
/*   fprintf(outfile, "MUZ_AOs\n"); */
94
 
/*   print_mat(MUZ_AO, nao, nao, outfile); */
95
 
 
96
 
  MUX_MO = block_matrix(nso,nso);
97
 
  MUY_MO = block_matrix(nso,nso);
98
 
  MUZ_MO = block_matrix(nso,nso);
99
 
 
100
 
  /*** Transform the AO dipole integrals to the SO basis ***/
101
 
  X = block_matrix(nso,nao); /* just a temporary matrix */
102
 
 
103
 
  C_DGEMM('n','n',nso,nao,nao,1,&(usotao[0][0]),nao,&(MUX_AO[0][0]),nao,
104
 
          0,&(X[0][0]),nao);
105
 
  C_DGEMM('n','t',nso,nso,nao,1,&(X[0][0]),nao,&(usotao[0][0]),nao,
106
 
          0,&(MUX_MO[0][0]),nso);
107
 
 
108
 
  C_DGEMM('n','n',nso,nao,nao,1,&(usotao[0][0]),nao,&(MUY_AO[0][0]),nao,
109
 
          0,&(X[0][0]),nao);
110
 
  C_DGEMM('n','t',nso,nso,nao,1,&(X[0][0]),nao,&(usotao[0][0]),nao,
111
 
          0,&(MUY_MO[0][0]),nso);
112
 
 
113
 
  C_DGEMM('n','n',nso,nao,nao,1,&(usotao[0][0]),nao,&(MUZ_AO[0][0]),nao,
114
 
          0,&(X[0][0]),nao);
115
 
  C_DGEMM('n','t',nso,nso,nao,1,&(X[0][0]),nao,&(usotao[0][0]),nao,
116
 
          0,&(MUZ_MO[0][0]),nso);
117
 
 
118
 
  free(mu_x_ints); 
119
 
  free(mu_y_ints); 
120
 
  free(mu_z_ints);
121
 
  free_block(X);
122
 
  free_block(usotao);
123
 
  free_block(MUX_AO);
124
 
  free_block(MUY_AO);
125
 
  free_block(MUZ_AO);
126
 
 
127
 
  /*** Transform the SO dipole integrals to the MO basis ***/
128
 
 
129
 
  X = block_matrix(nmo,nmo); /* just a temporary matrix */
130
 
 
131
 
  C_DGEMM('t','n',nmo,nmo,nmo,1,&(scf_qt[0][0]),nmo,&(MUX_MO[0][0]),nmo,
132
 
          0,&(X[0][0]),nmo);
133
 
  C_DGEMM('n','n',nmo,nmo,nmo,1,&(X[0][0]),nmo,&(scf_qt[0][0]),nmo,
134
 
          0,&(MUX_MO[0][0]),nmo);
135
 
 
136
 
  C_DGEMM('t','n',nmo,nmo,nmo,1,&(scf_qt[0][0]),nmo,&(MUY_MO[0][0]),nmo,
137
 
          0,&(X[0][0]),nmo);
138
 
  C_DGEMM('n','n',nmo,nmo,nmo,1,&(X[0][0]),nmo,&(scf_qt[0][0]),nmo,
139
 
          0,&(MUY_MO[0][0]),nmo);
140
 
 
141
 
  C_DGEMM('t','n',nmo,nmo,nmo,1,&(scf_qt[0][0]),nmo,&(MUZ_MO[0][0]),nmo,
142
 
          0,&(X[0][0]),nmo);
143
 
  C_DGEMM('n','n',nmo,nmo,nmo,1,&(X[0][0]),nmo,&(scf_qt[0][0]),nmo,
144
 
          0,&(MUZ_MO[0][0]),nmo);
145
 
 
146
 
  free_block(scf_qt);
147
 
  free_block(X);
148
 
 
149
 
  moinfo.dip = (double ***) malloc(3 * sizeof(double **));
150
 
  moinfo.dip[0] = block_matrix(nao, nao);
151
 
  moinfo.dip[1] = block_matrix(nao, nao);
152
 
  moinfo.dip[2] = block_matrix(nao, nao);
153
 
 
154
 
  for(i=0; i<nmo; i++) 
155
 
    for(j=0; j<nmo; j++) {
156
 
      moinfo.dip[0][i][j] = MUX_MO[i][j];
157
 
      moinfo.dip[1][i][j] = MUY_MO[i][j];
158
 
      moinfo.dip[2][i][j] = MUZ_MO[i][j];
159
 
    }
160
 
 
161
 
  free(ioff);
162
 
  free_block(MUX_MO); 
163
 
  free_block(MUY_MO); 
164
 
  free_block(MUZ_MO);
165
 
 
166
 
  return;
167
 
}