1 /**********************************************************************
2   Hamiltonian_Cluster.c:
3 
4      Hamiltonian_Cluster.c is a subroutine to make a Hamiltonian matrix
5      for cluster or molecular systems.
6 
7   Log of Hamiltonian_Cluster.c:
8 
9      22/Nov/2001  Released by T.Ozaki
10 
11 ***********************************************************************/
12 
13 #include <stdio.h>
14 #include <stdlib.h>
15 #include <math.h>
16 #include "openmx_common.h"
17 #include "mpi.h"
18 
19 
Hamiltonian_Cluster(double **** RH,double ** H,int * MP)20 void Hamiltonian_Cluster(double ****RH, double **H, int *MP)
21 {
22   int i,j,k;
23   int MA_AN,GA_AN,LB_AN,GB_AN,AN;
24   int wanA,wanB,tnoA,tnoB,Anum,Bnum,NUM;
25   int num,tnum,num_orbitals;
26   int ID,myid,numprocs,tag=999;
27   int *My_NZeros;
28   int *is1,*ie1,*is2;
29   int *My_Matomnum,*order_GA;
30   double *H1,sum;
31 
32   MPI_Status stat;
33   MPI_Request request;
34 
35   /* MPI */
36 
37   MPI_Comm_size(mpi_comm_level1,&numprocs);
38   MPI_Comm_rank(mpi_comm_level1,&myid);
39   MPI_Barrier(mpi_comm_level1);
40 
41   /* allocation of arrays */
42 
43   My_NZeros = (int*)malloc(sizeof(int)*numprocs);
44   My_Matomnum = (int*)malloc(sizeof(int)*numprocs);
45   is1 = (int*)malloc(sizeof(int)*numprocs);
46   ie1 = (int*)malloc(sizeof(int)*numprocs);
47   is2 = (int*)malloc(sizeof(int)*numprocs);
48   order_GA = (int*)malloc(sizeof(int)*(atomnum+2));
49 
50   /* find my total number of non-zero elements in myid */
51 
52   My_NZeros[myid] = 0;
53   for (MA_AN=1; MA_AN<=Matomnum; MA_AN++){
54     GA_AN = M2G[MA_AN];
55     wanA = WhatSpecies[GA_AN];
56     tnoA = Spe_Total_CNO[wanA];
57 
58     num = 0;
59     for (LB_AN=0; LB_AN<=FNAN[GA_AN]; LB_AN++){
60       GB_AN = natn[GA_AN][LB_AN];
61       wanB = WhatSpecies[GB_AN];
62       tnoB = Spe_Total_CNO[wanB];
63       num += tnoB;
64     }
65 
66     My_NZeros[myid] += tnoA*num;
67   }
68 
69   for (ID=0; ID<numprocs; ID++){
70     MPI_Bcast(&My_NZeros[ID],1,MPI_INT,ID,mpi_comm_level1);
71   }
72 
73   tnum = 0;
74   for (ID=0; ID<numprocs; ID++){
75     tnum += My_NZeros[ID];
76   }
77 
78   is1[0] = 0;
79   ie1[0] = My_NZeros[0] - 1;
80 
81   for (ID=1; ID<numprocs; ID++){
82     is1[ID] = ie1[ID-1] + 1;
83     ie1[ID] = is1[ID] + My_NZeros[ID] - 1;
84   }
85 
86   /* set is2 and order_GA */
87 
88   My_Matomnum[myid] = Matomnum;
89   for (ID=0; ID<numprocs; ID++){
90     MPI_Bcast(&My_Matomnum[ID],1,MPI_INT,ID,mpi_comm_level1);
91   }
92 
93   is2[0] = 1;
94   for (ID=1; ID<numprocs; ID++){
95     is2[ID] = is2[ID-1] + My_Matomnum[ID-1];
96   }
97 
98   for (MA_AN=1; MA_AN<=Matomnum; MA_AN++){
99     order_GA[is2[myid]+MA_AN-1] = M2G[MA_AN];
100   }
101 
102   for (ID=0; ID<numprocs; ID++){
103     MPI_Bcast(&order_GA[is2[ID]],My_Matomnum[ID],MPI_INT,ID,mpi_comm_level1);
104   }
105 
106   /* set MP */
107 
108   Anum = 1;
109   for (i=1; i<=atomnum; i++){
110     MP[i] = Anum;
111     wanA = WhatSpecies[i];
112     Anum += Spe_Total_CNO[wanA];
113   }
114   NUM = Anum - 1;
115 
116   /* set H1 */
117 
118   H1 = (double*)malloc(sizeof(double)*(tnum+1));
119 
120   k = is1[myid];
121   for (MA_AN=1; MA_AN<=Matomnum; MA_AN++){
122     GA_AN = M2G[MA_AN];
123     wanA = WhatSpecies[GA_AN];
124     tnoA = Spe_Total_CNO[wanA];
125     for (i=0; i<tnoA; i++){
126       for (LB_AN=0; LB_AN<=FNAN[GA_AN]; LB_AN++){
127         GB_AN = natn[GA_AN][LB_AN];
128         wanB = WhatSpecies[GB_AN];
129         tnoB = Spe_Total_CNO[wanB];
130         for (j=0; j<tnoB; j++){
131           H1[k] = RH[MA_AN][LB_AN][i][j];
132           k++;
133 	}
134       }
135     }
136   }
137 
138   /* MPI H1 */
139 
140   for (ID=0; ID<numprocs; ID++){
141     k = is1[ID];
142     MPI_Bcast(&H1[k], My_NZeros[ID], MPI_DOUBLE, ID, mpi_comm_level1);
143   }
144 
145   /* H1 -> H */
146 
147   H[0][0] = NUM;
148   for (i=1; i<=NUM; i++){
149     for (j=1; j<=NUM; j++){
150       H[i][j] = 0.0;
151     }
152   }
153 
154   k = 0;
155   for (AN=1; AN<=atomnum; AN++){
156     GA_AN = order_GA[AN];
157     wanA = WhatSpecies[GA_AN];
158     tnoA = Spe_Total_CNO[wanA];
159     Anum = MP[GA_AN];
160 
161     for (i=0; i<tnoA; i++){
162 
163       for (LB_AN=0; LB_AN<=FNAN[GA_AN]; LB_AN++){
164         GB_AN = natn[GA_AN][LB_AN];
165         wanB = WhatSpecies[GB_AN];
166         tnoB = Spe_Total_CNO[wanB];
167         Bnum = MP[GB_AN];
168 
169         for (j=0; j<tnoB; j++){
170           H[Anum+i][Bnum+j] += H1[k];
171           k++;
172 	}
173       }
174     }
175   }
176 
177   /* freeing of arrays */
178 
179   free(My_NZeros);
180   free(My_Matomnum);
181   free(is1);
182   free(ie1);
183   free(is2);
184   free(order_GA);
185   free(H1);
186 }
187