1 /**********************************************************************
2 Hamiltonian_Band_NC.c:
3
4 Hamiltonian_Band_NC.c is a subroutine to make a spin non-collinear
5 Hamiltonian matrix for a periodic boundary system using Bloch theorem.
6
7 Log of Hamiltonian_Band_NC.c:
8
9 24/Dec/2003 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
Hamiltonian_Band_NC(int Host_ID1,double ***** RH,double ***** IH,dcomplex ** H,int * MP,double k1,double k2,double k3)19 void Hamiltonian_Band_NC(int Host_ID1,
20 double *****RH, double *****IH,
21 dcomplex **H, int *MP,
22 double k1, double k2, double k3)
23 {
24 static int firsttime=1;
25 int i,j,k,wanA,wanB,tnoA,tnoB,Anum,Bnum;
26 int NUM,MA_AN,GA_AN,LB_AN,GB_AN;
27 int l1,l2,l3,Rn,n2;
28 double *tmp_array1,*tmp_array2;
29 double **H11r,**H11i;
30 double **H22r,**H22i;
31 double **H12r,**H12i;
32 double kRn,si,co,h;
33 int ID,myid,numprocs,tag=999;
34
35 MPI_Status stat;
36 MPI_Request request;
37
38 /* MPI */
39 MPI_Comm_size(mpi_comm_level1,&numprocs);
40 MPI_Comm_rank(mpi_comm_level1,&myid);
41 MPI_Barrier(mpi_comm_level1);
42
43 /* set MP */
44 Anum = 1;
45 for (i=1; i<=atomnum; i++){
46 MP[i] = Anum;
47 wanA = WhatSpecies[i];
48 tnoA = Spe_Total_CNO[wanA];
49 Anum = Anum + tnoA;
50 }
51 NUM = Anum - 1;
52 n2 = NUM + 2;
53
54 /*******************************************
55 allocation of H11r, H11i,
56 H22r, H22i,
57 H12r, H12i
58 *******************************************/
59
60 H11r = (double**)malloc(sizeof(double*)*n2);
61 for (i=0; i<n2; i++){
62 H11r[i] = (double*)malloc(sizeof(double)*n2);
63 }
64
65 H11i = (double**)malloc(sizeof(double*)*n2);
66 for (i=0; i<n2; i++){
67 H11i[i] = (double*)malloc(sizeof(double)*n2);
68 }
69
70 H22r = (double**)malloc(sizeof(double*)*n2);
71 for (i=0; i<n2; i++){
72 H22r[i] = (double*)malloc(sizeof(double)*n2);
73 }
74
75 H22i = (double**)malloc(sizeof(double*)*n2);
76 for (i=0; i<n2; i++){
77 H22i[i] = (double*)malloc(sizeof(double)*n2);
78 }
79
80 H12r = (double**)malloc(sizeof(double*)*n2);
81 for (i=0; i<n2; i++){
82 H12r[i] = (double*)malloc(sizeof(double)*n2);
83 }
84
85 H12i = (double**)malloc(sizeof(double*)*n2);
86 for (i=0; i<n2; i++){
87 H12i[i] = (double*)malloc(sizeof(double)*n2);
88 }
89
90 /****************************************************
91 PrintMemory
92 ****************************************************/
93
94 if (firsttime){
95 PrintMemory("Hamiltonian_Band: H11r",sizeof(double)*n2*n2,NULL);
96 PrintMemory("Hamiltonian_Band: H11i",sizeof(double)*n2*n2,NULL);
97 PrintMemory("Hamiltonian_Band: H22r",sizeof(double)*n2*n2,NULL);
98 PrintMemory("Hamiltonian_Band: H22i",sizeof(double)*n2*n2,NULL);
99 PrintMemory("Hamiltonian_Band: H12r",sizeof(double)*n2*n2,NULL);
100 PrintMemory("Hamiltonian_Band: H12i",sizeof(double)*n2*n2,NULL);
101 }
102
103 /* for PrintMemory */
104 firsttime=0;
105
106 /****************************************************
107 set Hamiltonian
108 ****************************************************/
109
110 H[0][0].r = 2.0*NUM;
111 for (i=1; i<=NUM; i++){
112 for (j=1; j<=NUM; j++){
113 H11r[i][j] = 0.0;
114 H11i[i][j] = 0.0;
115 H22r[i][j] = 0.0;
116 H22i[i][j] = 0.0;
117 H12r[i][j] = 0.0;
118 H12i[i][j] = 0.0;
119 }
120 }
121
122 for (MA_AN=1; MA_AN<=Matomnum; MA_AN++){
123 GA_AN = M2G[MA_AN];
124 wanA = WhatSpecies[GA_AN];
125 tnoA = Spe_Total_CNO[wanA];
126 Anum = MP[GA_AN];
127
128 for (LB_AN=0; LB_AN<=FNAN[GA_AN]; LB_AN++){
129 GB_AN = natn[GA_AN][LB_AN];
130 Rn = ncn[GA_AN][LB_AN];
131 wanB = WhatSpecies[GB_AN];
132 tnoB = Spe_Total_CNO[wanB];
133
134 /*
135 kRn = k1*( (double)atv_ijk[Rn][1] + Cell_Gxyz[GB_AN][1] - Cell_Gxyz[GA_AN][1] )
136 + k2*( (double)atv_ijk[Rn][2] + Cell_Gxyz[GB_AN][2] - Cell_Gxyz[GA_AN][2] )
137 + k3*( (double)atv_ijk[Rn][3] + Cell_Gxyz[GB_AN][3] - Cell_Gxyz[GA_AN][3] );
138 */
139
140 l1 = atv_ijk[Rn][1];
141 l2 = atv_ijk[Rn][2];
142 l3 = atv_ijk[Rn][3];
143 kRn = k1*(double)l1 + k2*(double)l2 + k3*(double)l3;
144
145 si = sin(2.0*PI*kRn);
146 co = cos(2.0*PI*kRn);
147 Bnum = MP[GB_AN];
148
149 /* non-spin-orbit coupling and non-LDA+U */
150 if (SO_switch==0 && Hub_U_switch==0 && Constraint_NCS_switch==0
151 && Zeeman_NCS_switch==0 && Zeeman_NCO_switch==0){
152
153 for (i=0; i<=(tnoA-1); i++){
154 for (j=0; j<=(tnoB-1); j++){
155 H11r[Anum+i][Bnum+j] += co*RH[0][MA_AN][LB_AN][i][j];
156 H11i[Anum+i][Bnum+j] += si*RH[0][MA_AN][LB_AN][i][j];
157 H22r[Anum+i][Bnum+j] += co*RH[1][MA_AN][LB_AN][i][j];
158 H22i[Anum+i][Bnum+j] += si*RH[1][MA_AN][LB_AN][i][j];
159 H12r[Anum+i][Bnum+j] += co*RH[2][MA_AN][LB_AN][i][j] - si*RH[3][MA_AN][LB_AN][i][j];
160 H12i[Anum+i][Bnum+j] += si*RH[2][MA_AN][LB_AN][i][j] + co*RH[3][MA_AN][LB_AN][i][j];
161 }
162 }
163 }
164
165 /* spin-orbit coupling or LDA+U */
166 else {
167 for (i=0; i<=(tnoA-1); i++){
168 for (j=0; j<=(tnoB-1); j++){
169 H11r[Anum+i][Bnum+j] += co*RH[0][MA_AN][LB_AN][i][j] - si*IH[0][MA_AN][LB_AN][i][j];
170 H11i[Anum+i][Bnum+j] += si*RH[0][MA_AN][LB_AN][i][j] + co*IH[0][MA_AN][LB_AN][i][j];
171 H22r[Anum+i][Bnum+j] += co*RH[1][MA_AN][LB_AN][i][j] - si*IH[1][MA_AN][LB_AN][i][j];
172 H22i[Anum+i][Bnum+j] += si*RH[1][MA_AN][LB_AN][i][j] + co*IH[1][MA_AN][LB_AN][i][j];
173 H12r[Anum+i][Bnum+j] += co*RH[2][MA_AN][LB_AN][i][j] - si*(RH[3][MA_AN][LB_AN][i][j]
174 + IH[2][MA_AN][LB_AN][i][j]);
175 H12i[Anum+i][Bnum+j] += si*RH[2][MA_AN][LB_AN][i][j] + co*(RH[3][MA_AN][LB_AN][i][j]
176 + IH[2][MA_AN][LB_AN][i][j]);
177 }
178 }
179 }
180
181 }
182 }
183
184 /******************************************************
185 MPI: H11r, H11i
186 H22r, H22i
187 H12r, H12i
188 ******************************************************/
189
190 tmp_array1 = (double*)malloc(sizeof(double)*6*n2*List_YOUSO[7]);
191 tmp_array2 = (double*)malloc(sizeof(double)*6*n2*List_YOUSO[7]);
192
193 if (myid!=Host_ID1){
194
195 for (MA_AN=1; MA_AN<=Matomnum; MA_AN++){
196 GA_AN = M2G[MA_AN];
197 wanA = WhatSpecies[GA_AN];
198 tnoA = Spe_Total_CNO[wanA];
199 Anum = MP[GA_AN];
200
201 k = 0;
202
203 for (i=0; i<tnoA; i++){
204 for (j=0; j<n2; j++){
205 tmp_array1[k] = H11r[Anum+i][j];
206 k++;
207 }
208 }
209
210 for (i=0; i<tnoA; i++){
211 for (j=0; j<n2; j++){
212 tmp_array1[k] = H11i[Anum+i][j];
213 k++;
214 }
215 }
216
217 for (i=0; i<tnoA; i++){
218 for (j=0; j<n2; j++){
219 tmp_array1[k] = H22r[Anum+i][j];
220 k++;
221 }
222 }
223
224 for (i=0; i<tnoA; i++){
225 for (j=0; j<n2; j++){
226 tmp_array1[k] = H22i[Anum+i][j];
227 k++;
228 }
229 }
230
231 for (i=0; i<tnoA; i++){
232 for (j=0; j<n2; j++){
233 tmp_array1[k] = H12r[Anum+i][j];
234 k++;
235 }
236 }
237
238 for (i=0; i<tnoA; i++){
239 for (j=0; j<n2; j++){
240 tmp_array1[k] = H12i[Anum+i][j];
241 k++;
242 }
243 }
244
245 tag = 999;
246 MPI_Isend(&tmp_array1[0], 6*tnoA*n2, MPI_DOUBLE, Host_ID1,
247 tag, mpi_comm_level1, &request);
248 MPI_Wait(&request,&stat);
249
250 }
251
252 }
253
254 else{
255
256 for (GA_AN=1; GA_AN<=atomnum; GA_AN++){
257 wanA = WhatSpecies[GA_AN];
258 tnoA = Spe_Total_CNO[wanA];
259 Anum = MP[GA_AN];
260 ID = G2ID[GA_AN];
261 if (ID!=Host_ID1){
262
263 tag = 999;
264 MPI_Recv(&tmp_array2[0], 6*tnoA*n2, MPI_DOUBLE, ID, tag, mpi_comm_level1, &stat);
265
266 k = 0;
267
268 for (i=0; i<tnoA; i++){
269 for (j=0; j<n2; j++){
270 H11r[Anum+i][j] = tmp_array2[k];
271 k++;
272 }
273 }
274
275 for (i=0; i<tnoA; i++){
276 for (j=0; j<n2; j++){
277 H11i[Anum+i][j] = tmp_array2[k];
278 k++;
279 }
280 }
281
282 for (i=0; i<tnoA; i++){
283 for (j=0; j<n2; j++){
284 H22r[Anum+i][j] = tmp_array2[k];
285 k++;
286 }
287 }
288
289 for (i=0; i<tnoA; i++){
290 for (j=0; j<n2; j++){
291 H22i[Anum+i][j] = tmp_array2[k];
292 k++;
293 }
294 }
295
296 for (i=0; i<tnoA; i++){
297 for (j=0; j<n2; j++){
298 H12r[Anum+i][j] = tmp_array2[k];
299 k++;
300 }
301 }
302
303 for (i=0; i<tnoA; i++){
304 for (j=0; j<n2; j++){
305 H12i[Anum+i][j] = tmp_array2[k];
306 k++;
307 }
308 }
309
310 }
311 }
312
313 }
314
315 free(tmp_array1);
316 free(tmp_array2);
317
318 /******************************************************
319 the full complex matrix of H
320 ******************************************************/
321
322 for (i=1; i<=NUM; i++){
323 for (j=1; j<=NUM; j++){
324 H[i ][j ].r = H11r[i][j];
325 H[i ][j ].i = H11i[i][j];
326 H[i+NUM][j+NUM].r = H22r[i][j];
327 H[i+NUM][j+NUM].i = H22i[i][j];
328 H[i ][j+NUM].r = H12r[i][j];
329 H[i ][j+NUM].i = H12i[i][j];
330 H[j+NUM][i ].r = H[i][j+NUM].r;
331 H[j+NUM][i ].i = -H[i][j+NUM].i;
332 }
333 }
334
335 /****************************************************
336 free arrays
337 ****************************************************/
338
339 for (i=0; i<n2; i++){
340 free(H11r[i]);
341 }
342 free(H11r);
343
344 for (i=0; i<n2; i++){
345 free(H11i[i]);
346 }
347 free(H11i);
348
349 for (i=0; i<n2; i++){
350 free(H22r[i]);
351 }
352 free(H22r);
353
354 for (i=0; i<n2; i++){
355 free(H22i[i]);
356 }
357 free(H22i);
358
359 for (i=0; i<n2; i++){
360 free(H12r[i]);
361 }
362 free(H12r);
363
364 for (i=0; i<n2; i++){
365 free(H12i[i]);
366 }
367 free(H12i);
368
369 }
370
371
372
373
374
375