Updated to the current version of MKL.
[ptslat.git] / input_pzo.f90
blob8a1b8cd4e1aaa5e44df444516d7fb689b8387889
1 MODULE Input_Data
2 IMPLICIT NONE
4 REAL, SAVE :: A_S,B_S,C_S ! -> Superlattice constants
5 REAL, DIMENSION(3), SAVE :: A1_S,A2_S,A3_S ! -> Superlattice vectors
6 INTEGER, SAVE :: NMin_X,NMax_X,NMin_Y,NMax_Y,NMin_Z,NMax_Z
7 REAL, SAVE :: X_Min,X_Max,Y_Min,Y_Max,Z_Min,Z_Max
8 REAL, SAVE :: X_Inc,Y_Inc,Z_Inc
9 INTEGER, SAVE :: XDim,YDim,ZDim
10 INTEGER, SAVE :: KCOOR
11 INTEGER,PRIVATE :: NDots_X,NDots_Y,NDots_Z
13 !!! Material Parameters
14 REAL, SAVE :: XLAMB,XMU,C13,C33,C11,EPSC,EPSA
15 !!! Calculation Constants
16 REAL, SAVE :: ETA1,ETA2,ETA_DIF,CN31_1,CN31_2,&
17 BISUM,BIZZ,BIAUX
19 !!! Normalized Dot Parameters
20 REAL, SAVE :: RC,ZC,D,RD,HD
22 !!! PIEZOELECTRIC CONSTANTS
23 REAL, PRIVATE :: PZE15,PZE31,PZE33,PSPB,PSPW,DIELEC
24 REAL, SAVE :: CI1,CI2,CI3
25 REAL, SAVE :: CNE2_1,CNE2_2,CNE_1,CNE_2,CSP,CPZ
26 REAL, SAVE :: CSPWL,CSPBR,CPZWL,CPZ_Q1,CPZ_Q2
28 CONTAINS
30 SUBROUTINE READ_INPUT()
32 Use Dot_Geometry
33 Use Auxiliar_Procedures, ONLY : AISO, MTYPE
35 IMPLICIT NONE
37 INTEGER :: G_Method
38 REAL :: D_Z
39 REAL,DIMENSION(3) :: V_AUX
43 !!!!! read data from data file !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
45 READ (*,*)
46 READ (*,*)
47 READ (*,*)
48 READ (*,*)
49 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
50 READ (*,*)
51 READ (*,*) DWL ! Wetting layer thickness (A)
52 READ (*,*) ISHAPE ! Shape of the dot
53 READ (*,*) HQD ! Height of the quantum dot (A)
54 READ (*,*) Rqd_Base ! Base Radius (A)
55 READ (*,*) Rqd_Top ! Top Radius (A)
56 !!! Superlattice Begins !!!!
57 READ (*,*)
58 READ (*,*) A_S,B_S,C_S
59 READ (*,*) A1_S(1),A1_S(2),A1_S(3)
60 READ (*,*) A2_S(1),A2_S(2),A2_S(3)
61 READ (*,*) A3_S(1),A3_S(2),A3_S(3)
62 READ (*,*) NDots_X,NDots_Y,NDots_Z
63 V_AUX=(/A_S,B_S,C_S/)
64 A1_S=A1_S*V_AUX
65 A2_S=A2_S*V_AUX
66 A3_S=A3_S*V_AUX
67 IF(MOD(NDots_X,2).EQ.0) THEN
68 NDots_X=NDots_X+1
69 WRITE(6,'(A,1X,I3)')"Warning: NDots_X Even -> NDots_X=",NDots_X
70 END IF
71 IF(MOD(NDots_Y,2).EQ.0) THEN
72 NDots_Y=NDots_Y+1
73 WRITE(6,'(A,1X,I3)')"Warning: NDots_Y Even -> NDots_Y=",NDots_Y
74 END IF
75 IF(MOD(NDots_Z,2).EQ.0) THEN
76 NDots_Z=NDots_Z+1
77 WRITE(6,'(A,1X,I3)')"Warning: NDots_Z Even -> NDots_Z=",NDots_Z
78 END IF
79 IF(NDots_X.EQ.1) THEN
80 NMin_X=0; NMax_X=0
81 ELSE
82 NMin_X=-(NDots_X-1)/2; NMax_X=-NMin_X
83 END IF
84 IF(NDots_Y.EQ.1) THEN
85 NMin_Y=0; NMax_Y=0
86 ELSE
87 NMin_Y=-(NDots_Y-1)/2; NMax_Y=-NMin_Y
88 END IF
89 IF(NDots_Z.EQ.1) THEN
90 NMin_Z=0; NMax_Z=0
91 ELSE
92 NMin_Z=-(NDots_Z-1)/2; NMax_Z=-NMin_Z
93 END IF
94 !!! Grid Begins !!!!
95 READ (*,*)
96 READ (*,*) G_Method
97 READ (*,*) X_Min,X_Max
98 READ (*,*) Y_Min,Y_Max
99 READ (*,*) Z_Min,Z_Max
100 READ (*,*) XDim,YDim,ZDim
102 IF(G_Method.EQ.1) THEN
103 IF(XDim.EQ.1) THEN
104 X_Max=X_Min
105 ELSE
106 X_Max=A1_S(1)/2.E0
107 X_Min=0.E0
108 END IF
109 IF(YDim.EQ.1) THEN
110 Y_Max=Y_Min
111 ELSE
112 Y_Max=A2_S(2)/2.E0
113 Y_Min=0.E0
114 END IF
115 IF(ZDim.EQ.1) THEN
116 Z_Max=Z_Min
117 ELSE
118 Z_Max=A3_S(3)/2.E0
119 D_Z=Z_Max-(Hqd+DWL)/2.E0
120 Z_Max=Hqd+D_Z
121 Z_Min=-DWL-D_Z
122 END IF
123 END IF
125 IF (XDim.EQ.1) THEN
126 X_Inc=0.E0
127 ELSE
128 X_Inc=(X_Max-X_Min)/REAL(XDim-1)
129 END IF
130 IF (YDim.EQ.1) THEN
131 Y_Inc=0.E0
132 ELSE
133 Y_Inc=(Y_Max-Y_Min)/REAL(YDim-1)
134 END IF
135 IF (ZDim.EQ.1) THEN
136 Z_Inc=0.E0
137 ELSE
138 Z_Inc=(Z_Max-Z_Min)/REAL(ZDim-1)
139 END IF
141 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
143 READ (*,*)
144 READ (*,*) C11,C13,C33 ! Elastic modulus
145 READ (*,*) XLAMB,XMU ! Lame constants
146 READ (*,*) EPSA ! Misfit strain: EPS0
147 READ (*,*) EPSC ! Misfit strain: EPSC
148 READ (*,*) PZE15,PZE31,PZE33 !Piezoelectric Constants
149 READ (*,*) PSPB,PSPW ! Spontaneous Polarization
150 READ (*,*) DIELEC ! Barrier Dielectric Constant
151 READ (*,*) MTYPE ! MTYPE: 1-> WZI, 2-> WZA
152 IF (MTYPE.LT.2) THEN
153 ! C13=XLAMB; C11=XLAMB+2.E0*XMU; C33=C11
154 ! EPSC=EPSA
155 AISO=1
156 ! WRITE(6,*)"ERROR: WZ ISO NOT IMPLEMENTED"
157 ! STOP
158 ELSE
159 AISO=0
160 END IF
161 READ (*,*) KCOOR ! KCOOR=0 Cartesian Coordinates
162 ! KCOOR=1 Cylindrical Coordinates
163 IF(KCOOR.EQ.1.AND.YDIM.GT.1) THEN
164 WRITE(6,*)"ERROR: Selected coordinates are cylindrical and YDIM.GT.1"
165 WRITE(6,*)" Exiting Program"
166 STOP
167 END IF
168 READ (*,*)
169 !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
172 RETURN
174 END SUBROUTINE READ_INPUT
176 SUBROUTINE CONSTANTS( )
177 Use Dot_Geometry, ONLY: DWL
178 Use Auxiliar_Procedures, ONLY : AISO
179 IMPLICIT NONE
181 REAL :: ALPHA,BETA,GAMA
182 REAL :: R,S,A1,A2,B1,B2
183 REAL :: CN3_1,CN3_2,P1,P2,P3
184 REAL :: ETA_AUX, PI
186 PI=4.E0*ATAN(1.E0)
188 !!!! Remember: XLAMB=C12, XMU=C44 !!!!!
190 !!! ISO CONSTANTS
191 BISUM = 2.*EPSA ! Biaxial strain components
192 BIZZ = -2.*EPSA*XLAMB/(XLAMB+2.*XMU)
193 BIAUX = EPSA-BIZZ
195 R=(C11+XLAMB+C13)*EPSA+C13*(EPSC-EPSA)
196 S=(C33+C13-XLAMB-C11)*EPSA+(C33-C13)*(EPSC-EPSA)
198 ALPHA=-XMU*C33
200 R=R/ALPHA; S=S/ALPHA
202 A1=R*(C33-XMU-C13)-S*(C13+XMU)
203 A2=R*XMU
205 B1=S*(C13+2.E0*XMU)-R*(C33-C13-2.E0*XMU)
206 B2=S*C11-R*(C13+2.E0*XMU-C11)
208 BETA=( C13*(C13+2.E0*XMU)-C11*C33 )/(-XMU*C33)
209 GAMA=C11/C33
211 CN3_1=A1+B1; CN3_2=A2+B2
213 P1=2.E0*PZE15; P2=PZE31; P3=PZE33
215 CI1=2.E0*EPSA*P2+EPSC*P3
216 CI2=P3
217 CI3=1.E0
219 CNE_1=(P1+P2-P3)*A1+(P1/2-P3)*B1
220 CNE_2=(P1+P2-P3)*A2+(P1/2-P3)*B2
222 ETA_AUX=SQRT(BETA**2 - 4.E0*GAMA)
223 ETA1=(-BETA+ETA_AUX)/2.E0
224 ETA2=(-BETA-ETA_AUX)/2.E0
225 IF(ETA1.GT.0.E0.OR.ETA2.GT.0.E0) THEN
226 WRITE(6,*)"ERROR: ETA1 AND ETA2 ARE GREATHER THAN ZERO"
227 STOP
228 END IF
229 ETA1=ABS(ETA1)
230 ETA2=ABS(ETA2)
232 ETA_DIF=ETA1-ETA2
234 CN31_1=(CN3_2/SQRT(ETA2)-CN3_1*SQRT(ETA2))
235 CN31_2=(CN3_2/SQRT(ETA1)-CN3_1*SQRT(ETA1))
237 CNE2_1=(1.E0-ETA1)*(CNE_2/SQRT(ETA2)-CNE_1*SQRT(ETA2))
238 CNE2_2=(1.E0-ETA2)*(CNE_2/SQRT(ETA1)-CNE_1*SQRT(ETA1))
240 CSP=(PSPW-PSPB)
241 !!! V=q/(4*PI*E0)*PZ=0.9*PZ eV
242 CPZ=0.9*(4.E0*PI/DIELEC)
243 CPZ_Q1=(BISUM*PZE31+BIZZ*PZE33)
244 CPZ_Q2=(BIAUX*(PZE33-(PZE31+2.E0*PZE15)))
246 !!! Wetting Layer Parameters
248 IF(NDots_z.GT.1) THEN
249 CSPWL=CSP*(C_S-DWL)/C_S
250 CSPBR=-CSP*(DWL/C_S)
251 ELSE
252 CSPWL=CSP
253 CSPBR=0.E0
254 END IF
256 IF(AISO.EQ.0) THEN
257 CPZWL=(PZE31-PZE33*(C13/C33))*BISUM
258 ELSE
259 CPZWL=(BISUM*PZE31+BIZZ*PZE33)
260 END IF
262 RETURN
264 END SUBROUTINE CONSTANTS
267 END MODULE INPUT_DATA