4
#include <libciomr/libciomr.h>
11
void build_so_classes()
16
int class,lmax,irr,ao,symop,char_i;
17
int ao_first, ao_last, ao_max, ao_irr;
18
int class_first, class_last, uc;
19
int shell, shell_first, shell_last;
20
int **class_num_angmom;
22
int so_cnt, so_blk_lim, so_blk_off, cnt;
24
int basisfn_num_in_symblk;
25
int *first_basisfn_in_symblk_offset;
26
int *first_basisfn_in_symblk;
29
/*----------------------------------------------------------------
30
Compute maximum angular momentum number for each class of atoms
31
----------------------------------------------------------------*/
33
max_angmom_class = init_int_array(num_classes);
34
class_num_angmom = init_int_matrix(num_classes,max_angmom+1);
35
for(i=0;i<num_atoms;i++) {
36
kfirst = first_shell_on_atom[i];
37
klast = kfirst + nshells_per_atom[i];
38
class = atom_class[i];
39
for(j=kfirst;j<klast;j++) {
40
max_angmom_class[class] = (max_angmom_class[class] > shell_ang_mom[j]) ? max_angmom_class[class] : shell_ang_mom[j];
41
class_num_angmom[class][shell_ang_mom[j]]++;
45
/*-----------------------
46
Allocate global arrays
47
-----------------------*/
49
num_cart_so_in_class = (int ***) malloc(num_classes*sizeof(int **));
50
for(class=0;class<num_classes;class++)
51
num_cart_so_in_class[class] = init_int_matrix(max_angmom_class[class]+1,nirreps);
52
class_so_coeff = (double ****) malloc(num_classes*sizeof(double ***));
53
for(i=0;i<num_classes;i++)
54
class_so_coeff[i] = (double ***) malloc((max_angmom_class[i]+1)*sizeof(double **));
56
num_pureang_so_in_class = (int ***) malloc(num_classes*sizeof(int **));
57
for(class=0;class<num_classes;class++)
58
num_pureang_so_in_class[class] = init_int_matrix(max_angmom_class[class]+1,nirreps);
62
/*----------------------
64
-----------------------*/
66
coeff_irr = init_int_matrix(nirreps,ioff[max_angmom+1]);
69
/*------------------------------------------------------------------------------------------------------------------------
70
Loop over each unique class, each angular momentum number within class, irreps, and basis function types
71
and 1) count number of SO in each representation num_cart_so_in_class[class][l][irr] generated by a shell of ang momentum l
72
on an atom of a unique class "uc" or its symmetry equivalent classes;
73
2) form matrices of coefficients and labels for classes equivalent to the unique class "uc"
74
------------------------------------------------------------------------------------------------------------------------*/
75
for(uc=0;uc<num_unique_classes;uc++) {
78
class_last = class + unique_class_degen[uc];
79
lmax = max_angmom_class[class];
80
for(l=0;l<=lmax;l++) {
82
for(symop=0;symop<nirreps;symop++)
83
if (class_orbit[class][symop] == class)
84
for(irr=0;irr<nirreps;irr++)
85
for(ao=0;ao<ao_max;ao++)
86
coeff_irr[irr][ao] += irr_char[irr][symop]*ao_type_transmat[l][symop][ao];
87
for(irr=0;irr<nirreps;irr++)
88
for(ao=0;ao<ao_max;ao++)
89
if (coeff_irr[irr][ao] > 0)
90
for(i=class_first;i<class_last;i++)
91
num_cart_so_in_class[i][l][irr] += 1;
94
for(irr=0;irr<nirreps;irr++)
95
for(i=class_first;i<class_last;i++)
96
num_pureang_so_in_class[i][l][irr] = num_cart_so_in_class[i][l][irr];
97
for(ao_irr=0;ao_irr<nirreps;ao_irr++)
98
if (num_cart_so[l][ao_irr] != 0)
99
for(ao=0;ao<ao_max;ao++)
100
if (ao_type_irr[l][ao] == ao_irr) {
101
for(irr=0;irr<nirreps;irr++)
102
if (coeff_irr[irr][ao] > 0)
103
for(i=class_first;i<class_last;i++)
104
num_pureang_so_in_class[i][l][irr] -= num_redun_so[l][ao_irr];
109
/*---------------------------------
110
Remove after testing is complete
111
---------------------------------*/
112
if (print_lvl >= DEBUGPRINT) {
113
fprintf(outfile," Class = %d l = %d\n Number of cartesian SOs in symmetry blocks:",class+1,l);
114
for(irr=0;irr<nirreps;irr++)
115
fprintf(outfile," %d",num_cart_so_in_class[class][l][irr]);
116
fprintf(outfile,"\n\n");
118
fprintf(outfile," Class = %d l = %d\n Number of pure angular momentum SOs in symmetry blocks:",class+1,l);
119
for(irr=0;irr<nirreps;irr++)
120
fprintf(outfile," %d",num_pureang_so_in_class[class][l][irr]);
121
fprintf(outfile,"\n\n");
125
/*Compute SO coefficients and labels for ALL classes - children of a unique class "class"*/
126
for(i=class_first;i<class_last;i++) {
127
class_so_coeff[i][l] = init_matrix(ao_max*unique_class_degen[uc],ao_max);
128
for(irr=0;irr<nirreps;irr++)
129
bzero((char *)coeff_irr[irr],sizeof(int)*ao_max);
130
for(symop=0;symop<nirreps;symop++)
131
if (class_orbit[class_first][symop] == i)
132
for(irr=0;irr<nirreps;irr++)
133
if (num_cart_so_in_class[class][l][irr] != 0) /*There are SOs in this class and ang. momentum type*/
134
for(ao=0;ao<ao_max;ao++)
135
coeff_irr[irr][ao] += irr_char[irr][symop]*ao_type_transmat[l][symop][ao];
137
for(irr=0;irr<nirreps;irr++)
138
if (num_cart_so_in_class[class][l][irr] != 0)
139
for(ao=0;ao<ao_max;ao++)
140
if (coeff_irr[irr][ao] != 0)
141
class_so_coeff[i][l][so_cnt++][ao] = sign(coeff_irr[irr][ao]);
143
for(irr=0;irr<nirreps;irr++)
144
bzero((char *)coeff_irr[irr],sizeof(int)*ao_max);
148
/*---------------------------------
149
Remove after testing is complete
150
---------------------------------*/
151
if (print_lvl >= DEBUGPRINT)
152
for(uc=0;uc<num_unique_classes;uc++) {
155
class_last = class + unique_class_degen[uc];
156
lmax = max_angmom_class[class];
157
for(i=class_first;i<class_last;i++) {
158
fprintf(outfile,"\n Class = %d\n",i+1);
159
for(l=0;l<=lmax;l++) {
161
fprintf(outfile," l = %d\n",l);
162
for(j=0;j<ao_max*unique_class_degen[uc];j++) {
163
fprintf(outfile," ");
164
for(ao=0;ao<ao_max;ao++)
165
fprintf(outfile," %4.1lf",class_so_coeff[i][l][j][ao]);
166
fprintf(outfile,"\n");
added
removed
src/bin/input/build_so_classes.c
