File: TRAN_Calc_GC_LorR.c

package info (click to toggle)
openmx 3.2.4.dfsg-3
  • links: PTS, VCS
  • area: main
  • in suites: lenny, squeeze
  • size: 62,572 kB
  • ctags: 2,684
  • sloc: ansic: 130,666; python: 876; makefile: 560; xml: 63; perl: 18; sh: 4
file content (162 lines) | stat: -rw-r--r-- 4,216 bytes parent folder | download
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
/**********************************************************************
  TRNA_Calc_GC_LorR.c:
 
    TRAN_Calc_GC_LorR.c is a subroutine to calculate G_{C_L} or G_{C_R}
    in the non-equilibrium case.  

  Log of TRAN_Calc_GC_LorR.c:

     17/Dec/2005  Released by T.Ozaki, Taisuke Ozaki Copyright (C) 

***********************************************************************/

#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <math.h>

#ifdef nompi
#include "mimic_mpi.h"
#else
#include <mpi.h>
#endif

#include "tran_prototypes.h"
#include "lapack_prototypes.h"


/*
 *  calculate GCLorR
 *     GCLorR (w) = GC(w) ( w^* SC - HC - \Sigma_LorR(w^*) ) GC(w)^* 
 * 
 *
 *  no implicit variables 
 */


void TRAN_Calc_GC_LorR(
		       int iw_method,       /*  input */
                       dcomplex w,          /*  input */
                       double ChemP_e[2],   /*  input */
                       int nc,              /*  input */
                       int ne[2],           /*  input */
                       dcomplex *SigmaLorR, /*  input */
                       dcomplex *GC,        /*  input */
                       dcomplex *HCC,       /*  input */ 
                       dcomplex *SCC,       /*  input */ 
                       dcomplex *v1,        /* work, nc*nc */
                       dcomplex *GCLorR     /*  output */
                      )
#define GC_ref(i,j)  GC[nc*((j)-1)+(i)-1]
#define SigmaLorR_ref(i,j) SigmaLorR[nc*((j)-1)+(i)-1]
#define v1_ref(i,j)  v1[nc*((j)-1)+(i)-1] 
#define GCLorR_ref(i,j)  GCLorR[nc*((j)-1)+(i)-1]
#define HCC_ref(i,j) HCC[nc*((j)-1)+(i)-1]
#define SCC_ref(i,j) SCC[nc*((j)-1)+(i)-1]

{
  int i,j;
  int full_flag;
  int side;
  dcomplex alpha,beta;

  alpha.r = 1.0;
  alpha.i = 0.0;
  beta.r  = 0.0;
  beta.i  = 0.0;

  /* find a flag for how the H and S are divided */

  if (iw_method==3){
    if (ChemP_e[0]<ChemP_e[1])
      full_flag =  1;
    else 
      full_flag =  0;
  }
  else if (iw_method==4){
    if (ChemP_e[1]<ChemP_e[0])
      full_flag =  1;
    else 
      full_flag =  0;
  }

  if ( fabs(ChemP_e[1]-ChemP_e[0])<1.0e-50 ) full_flag = 2;

  /* GCLorR = -(\Sigma_LorR(w))^* */

  for (j=1; j<=nc; j++) {
    for (i=1; i<=nc; i++) {
      GCLorR_ref(i,j).r = -SigmaLorR_ref(i,j).r;
      GCLorR_ref(i,j).i =  SigmaLorR_ref(i,j).i;
    }
  }

  /* in case of the left side */

  if (iw_method==3){
    for (j=1; j<=(nc-ne[1]); j++) {
      for (i=1; i<=(nc-ne[1]); i++) {
        GCLorR_ref(i,j).r +=  w.r*SCC_ref(i,j).r - w.i*SCC_ref(i,j).i - HCC_ref(i,j).r; 
        GCLorR_ref(i,j).i += -w.i*SCC_ref(i,j).r - w.r*SCC_ref(i,j).i + HCC_ref(i,j).i; 
      }
    }
  }

  /* in case of the right side */

  else if (iw_method==4){
    for (j=(ne[0]+1); j<=nc; j++) {
      for (i=(ne[0]+1); i<=nc; i++) {
        GCLorR_ref(i,j).r +=  w.r*SCC_ref(i,j).r - w.i*SCC_ref(i,j).i - HCC_ref(i,j).r;
        GCLorR_ref(i,j).i += -w.i*SCC_ref(i,j).r - w.r*SCC_ref(i,j).i + HCC_ref(i,j).i;
      }
    }
  }

  /* correction of the central region */

  if (full_flag==0){
    for (j=(ne[0]+1); j<=(nc-ne[1]); j++) {
      for (i=(ne[0]+1); i<=(nc-ne[1]); i++) {
	GCLorR_ref(i,j).r -=  w.r*SCC_ref(i,j).r - w.i*SCC_ref(i,j).i - HCC_ref(i,j).r;
	GCLorR_ref(i,j).i -= -w.i*SCC_ref(i,j).r - w.r*SCC_ref(i,j).i + HCC_ref(i,j).i; 
      }
    }
  } 

  else if (full_flag==2){
    for (j=(ne[0]+1); j<=(nc-ne[1]); j++) {
      for (i=(ne[0]+1); i<=(nc-ne[1]); i++) {
	GCLorR_ref(i,j).r -= 0.5*(  w.r*SCC_ref(i,j).r - w.i*SCC_ref(i,j).i - HCC_ref(i,j).r );
	GCLorR_ref(i,j).i -= 0.5*( -w.i*SCC_ref(i,j).r - w.r*SCC_ref(i,j).i + HCC_ref(i,j).i );
      }
    }
  } 

  /* GC(w) ( w^* SC - HC -2 \Sigma_LorR(w^*) ) */

  F77_NAME(zgemm,ZGEMM)("N","N",&nc,&nc,&nc,&alpha, GC, &nc, GCLorR, &nc, &beta, v1, &nc);

  /* do complex conjugate */

  for (j=1;j<=nc;j++) {
    for (i=1;i<=nc;i++) {
      GC_ref(i,j).i = -GC_ref(i,j).i;
    }
  }

  /* GC(w) ( w^* SC - HC -2 \Sigma_LorR(w^*) ) GC(w)^* */
  
  F77_NAME(zgemm,ZGEMM)("N","N",&nc,&nc,&nc,&alpha, v1, &nc, GC, &nc, &beta, GCLorR, &nc);
  
  /* again do complex conjugate */
  
  for (j=1;j<=nc;j++) {
    for (i=1;i<=nc;i++) {
      GC_ref(i,j).i = -GC_ref(i,j).i;
    }
  }
}