← Back to branch summary

~ubuntu-branches/ubuntu/precise/psicode/precise

« back to all changes in this revision

Viewing changes to src/bin/stable/diag_A_UHF.cc

  • 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/stable/diag_A_UHF.cc
 
1
/*! \file
 
2
    \ingroup STABLE
 
3
    \brief Enter brief description of file here 
 
4
*/
 
5
#include <cstdio>
 
6
#include <cstdlib>
 
7
#include <libciomr/libciomr.h>
 
8
#include <libdpd/dpd.h>
 
9
#include <libqt/qt.h>
 
10
#include <psifiles.h>
 
11
#include "MOInfo.h"
 
12
#include "Params.h"
 
13
#define EXTERN
 
14
#include "globals.h"
 
15
 
 
16
namespace psi { namespace stable {
 
17
 
 
18
extern void follow_evec_UHF(double *vf, int dim_A, int dim_B);
 
19
extern void follow_evec_UHF2(double *vf, int dim_A, int dim_B);
 
20
 
 
21
void diag_A_UHF(void)
 
22
{
 
23
  int h, i, dim, dim_A, dim_B;
 
24
  int nroot, root;
 
25
  int ai, bj, ck, c, k, C, K, Ksym;
 
26
  double **A, *eps, **v, *evec_to_follow;
 
27
  char lbl[32];
 
28
  dpdbuf4 A_AA, A_BB, A_AB;
 
29
  dpdfile2 B;
 
30
 
 
31
  dpd_buf4_init(&A_AA, PSIF_MO_HESS, 0, 21, 21, 21, 21, 0, "A(AI,BJ)");
 
32
  dpd_buf4_init(&A_BB, PSIF_MO_HESS, 0, 31, 31, 31, 31, 0, "A(ai,bj)");
 
33
  dpd_buf4_init(&A_AB, PSIF_MO_HESS, 0, 21, 31, 21, 31, 0, "A(AI,bj)");
 
34
  for(h=0; h < moinfo.nirreps; h++) {
 
35
    dim_A = A_AA.params->rowtot[h];
 
36
    dim_B = A_BB.params->rowtot[h];
 
37
 
 
38
    dim = dim_A + dim_B;
 
39
    moinfo.rank[h] = dim;
 
40
 
 
41
    A = block_matrix(dim, dim);
 
42
    eps = init_array(dim);
 
43
    v = block_matrix(dim, dim);
 
44
 
 
45
    dpd_buf4_mat_irrep_init(&A_AA, h);
 
46
    dpd_buf4_mat_irrep_rd(&A_AA, h);
 
47
    for(ai=0; ai < dim_A; ai++)
 
48
      for(bj=0; bj < dim_A; bj++)
 
49
        A[ai][bj] = A_AA.matrix[h][ai][bj];
 
50
    dpd_buf4_mat_irrep_close(&A_AA, h);
 
51
 
 
52
    dpd_buf4_mat_irrep_init(&A_BB, h);
 
53
    dpd_buf4_mat_irrep_rd(&A_BB, h);
 
54
    for(ai=0; ai < dim_B; ai++)
 
55
      for(bj=0; bj < dim_B; bj++)
 
56
        A[ai+dim_A][bj+dim_A] = A_BB.matrix[h][ai][bj];
 
57
    dpd_buf4_mat_irrep_close(&A_BB, h);
 
58
 
 
59
    dpd_buf4_mat_irrep_init(&A_AB, h);
 
60
    dpd_buf4_mat_irrep_rd(&A_AB, h);
 
61
    for(ai=0; ai < dim_A; ai++)
 
62
      for(bj=0; bj < dim_B; bj++)
 
63
        A[ai][bj+dim_A] = A[bj+dim_A][ai] = A_AB.matrix[h][ai][bj];
 
64
    dpd_buf4_mat_irrep_close(&A_AB, h);
 
65
 
 
66
    sq_rsp(dim, dim, A, eps, 1, v, 1e-12);
 
67
 
 
68
    for(i=0; i < MIN0(moinfo.rank[h],5); i++)
 
69
      moinfo.A_evals[h][i] = eps[i];
 
70
 
 
71
    if (params.num_evecs_print > 0) {
 
72
      fprintf(outfile, "First %d eigenvectors for %s block:\n",
 
73
        MIN0(moinfo.rank[h], params.num_evecs_print), moinfo.labels[h]);
 
74
      eivout(v,eps,dim,MIN0(moinfo.rank[h],params.num_evecs_print),outfile); 
 
75
    }
 
76
 
 
77
    /* try to follow eigenvector downhill */
 
78
    if (h==0 && params.follow_instab) {
 
79
      if (eps[0] >= 0.0) 
 
80
        fprintf(outfile, "\nNo negative eigenvalues to follow.\n");
 
81
      else {
 
82
        evec_to_follow = init_array(dim);
 
83
        for (i=0; i<dim; i++) evec_to_follow[i] = v[i][0];
 
84
 
 
85
        fprintf(outfile, "\nAttempting to follow eigenvector using ");
 
86
        if (params.rotation_method == 0) {
 
87
          fprintf(outfile, "orbital rotation method.\n");
 
88
          follow_evec_UHF(evec_to_follow, dim_A, dim_B);
 
89
        }  
 
90
        else {
 
91
          fprintf(outfile, "antisymmetric matrix method.\n");
 
92
          follow_evec_UHF2(evec_to_follow, dim_A, dim_B);
 
93
        } 
 
94
        free(evec_to_follow);
 
95
      }
 
96
    }
 
97
 
 
98
    free(eps);
 
99
    free_block(v);
 
100
    free_block(A);
 
101
  }
 
102
 
 
103
  dpd_buf4_close(&A_AA);
 
104
  dpd_buf4_close(&A_BB);
 
105
  dpd_buf4_close(&A_AB);
 
106
 
 
107
}
 
108
 
 
109
 
 
110
}} // namespace psi::stable