#include "dims.h" SUBROUTINE LAMDAS use cons_module,only: dipmin,len1,len2,kmax,kmaxm1,kmaxp1, | imax,imaxp4,expz,boltz,avo,grav,gask,p0,rmassinv_o2, | rmassinv_o,rmassinv_n2 implicit none C **** C **** COMPUTE ION DRAG COEFFICIENTS C **** C **** MAJOR SPECIES AND IONS ARE NUMBERED THUS C **** O2 - 1 O2+ - 1 C **** O - 2 O+ - 2 C **** N2 - 3 N2+ - 3 C **** #include "params.h" #include "fcom.h" #include "index.h" #include "buff.h" #include "phys.h" #include "dynamo.h" #include "magfld.h" real :: rlatm,rlonm,dipmag,decmag,sndec,csdec,sn2dec,sncsdc, | dumdum COMMON/TRGM/RLATM(ZIMXP,ZJMX),RLONM(ZIMXP,ZJMX), 1 DIPMAG(ZIMXP,ZJMX),DECMAG(ZIMXP,ZJMX),SNDEC(ZIMXP,ZJMX), 2 CSDEC(ZIMXP,ZJMX),SN2DEC(ZIMXP,ZJMX),SNCSDC(ZIMXP,ZJMX), 3 DUMDUM(ZJMX,3) real :: rm1,rm2,rm3,v1,v2,v3,v4,v5,v6,t1,t2,t3,t4,t5,t6,t7 COMMON/VSCR/RM1(ZIMXP,ZKMXP,3,1),RM2(ZIMXP,ZKMXP,3,1), 1RM3(ZIMXP,ZKMXP,3,1),V1(ZIMXP,ZKMXP,1),V2(ZIMXP,ZKMXP,1), 2V3(ZIMXP,ZKMXP,1),V4(ZIMXP,ZKMXP,1),V5(ZIMXP,ZKMXP,1), 3V6(ZIMXP,ZKMXP,1),T1(ZIMXP),T2(ZIMXP),T3(ZIMXP),T4(ZIMXP), 4T5(ZIMXP),T6(ZIMXP),T7(ZIMXP) !$OMP THREADPRIVATE (/vscr/) !DIR$ TASKCOMMON vscr ! ! Local: real :: RMS(3) = (/32.,16.,30./) real,parameter :: QE=1.602E-19 integer :: l,m,i,k integer :: nps1k,nps2k,ntik,ntk,no2pk,nopk,nnopk,nmsk,ntek, | nlxyk,nzk,nuk,nvk,nwk,nlxxk,nlyyk ! C **** C **** T1 = BGAUSS C **** T2 = 9.6489E3 * BGAUS ( FOR OMEGA(I) ) C **** T3 = 1.7588028E7 * BGAUSS ( FOR OMEGA(E) ) C **** T4 = SIN(DIP) C **** DO 1 I=1,LEN1 T1(I)=BMOD(I,J) T2(I)=9.6489E3*T1(I) T3(I)=1.7588028E7*T1(I) ! ! DIPMIN is a parameter in init_mod.f (s.a., settei.F): if (abs(DIPMAG(I,J)) >= dipmin) then t4(i) = sin(dipmag(i,j)) else if(dipmag(i,j) >= 0.) then t4(i) = sin(dipmin) else t4(i) = sin(-dipmin) endif endif T5(I) = T4(I)**2 1 CONTINUE C **** C **** V1 = PSI(O2)/M(O2), V2 = PSI(O)/M(O), V1 = PSI2)/M(N2) C **** AT (K+1/2) C **** NPS1K=NJ+NPS NPS2K=NJ+NPS2 NTIK=NJ+NTI NTK=NJ+NT DO 2 I=1,LEN2 V1(I,1,1)=F(I,NPS1K)*rmassinv_o2 V2(I,1,1)=F(I,NPS2K)*rmassinv_o V3(I,1,1)=(1.-F(I,NPS1K)-F(I,NPS2K))*rmassinv_n2 C **** C **** V6 = TR = 0.5*(TI+TN) ( K+1/2 ) C **** V6(I,1,1)=0.5*(F(I,NTIK)+F(I,NTK)) 2 CONTINUE C **** C **** RM1 = NEW(L,M) C **** C **** GENERATE NUMERICAL FACTORS CALL NEW(V6,RM1) C **** MULTIPLY BY MAJOR SPECIES NUMBER DENSITIES DO 3 L=1,3 DO 3 M=1,3 DO 3 I=1,LEN2 RM1(I,1,L,M)=RM1(I,1,L,M)*V1(I,1,M) 3 CONTINUE C **** C **** V1(L) = SUM(M=1,3) RM1(L,M) ( K+1/2 ) C **** DO 4 L=1,3 DO 5 I=1,LEN2 V1(I,1,L)=0. 5 CONTINUE DO 4 M=1,3 DO 4 I=1,LEN2 V1(I,1,L)=V1(I,1,L)+RM1(I,1,L,M) 4 CONTINUE C **** C **** RM1(1) = N(O2+), RM1(2) = N(O+), RM1(1) = N(NO+) C **** ( AT (K+1/2) ) C **** NO2PK=NJ+NO2P NOPK=NJ+NOP NNOPK=NNOP DO 6 I=1,LEN2 RM1(I,1,1,1)=F(I,NO2PK) RM1(I,1,2,1)=F(I,NOPK) RM1(I,1,3,1)=F(I,NNOPK) 6 CONTINUE C **** C **** V6 = N * MBAR ( K+1/2 ) C **** NTK=NJ+NT-1 NMSK=NJ+NMS-1 DO 7 K=1,KMAX NTK=NTK+1 NMSK=NMSK+1 DO 7 I=1,LEN1 V6(I,K,1)=p0*expz(K)*.5*(F(I,NMSK)+F(I,NMSK+1))/(boltz*F(I,NTK)) 7 CONTINUE C **** C **** V1(L) = V1(L)*(N*MBAR)/OMEGAI(L) ( L=1,3 ) ( K=1/2 ) C **** DO 8 L=1,3 DO 8 K=1,KMAX DO 8 I=1,LEN1 V1(I,K,L)=V1(I,K,L)*V6(I,K,1)*RMS(L)/T2(I) 8 CONTINUE C **** C **** V4 = NU/OMEGA (E) C **** NPS1K=NJ+NPS NPS2K=NJ+NPS2 NTEK=NJ+NTE DO 9 I=1,LEN2 V5(I,1,1)=SQRT(F(I,NTEK)) V4(I,1,1)=V6(I,1,1)*(2.33E-11*(1.-F(I,NPS1K)-F(I,NPS2K))* | rmassinv_n2*F(I,NTEK)*(1.-1.21E-4*F(I,NTEK))+1.82E-10*F(I,NPS1K)* | rmassinv_o2*V5(I,1,1)*(1.+3.6E-2*V5(I,1,1))+8.9E-11*F(I,NPS2K)* | rmassinv_o*V5(I,1,1)*(1.+5.7E-4*F(I,NTEK))) 9 CONTINUE C **** C **** V4 = NU(EN)/OMEGAE C **** DO 10 K=1,KMAX DO 10 I=1,LEN1 V4(I,K,1)=V4(I,K,1)/T3(I) 10 CONTINUE C **** C **** RM2 = PEDERSEN FACTORS C **** RM3 = HALL FACTORS C **** DO 11 L=1,3 DO 11 I=1,LEN2 RM2(I,1,L,1)=V4(I,1,1)/(1.+V4(I,1,1)**2)+V1(I,1,L)/(1.+V1(I,1,L) 1**2) RM3(I,1,L,1)=1./(1.+V4(I,1,1)**2)-1./(1.+V1(I,1,L)**2) 11 CONTINUE C **** C **** V1 = SUM(L=1,3) PEDERSEN FACTORS * ION DENSITIES C **** V2 = SUM(L=1,3) HALL FACTORS * ION DENSITIES C **** DO 12 I=1,LEN2 V1(I,1,1)=0. V2(I,1,1)=0. 12 CONTINUE DO 13 L=1,3 DO 13 I=1,LEN2 V1(I,1,1)=V1(I,1,1)+RM2(I,1,L,1)*RM1(I,1,L,1) V2(I,1,1)=V2(I,1,1)+RM3(I,1,L,1)*RM1(I,1,L,1) 13 CONTINUE C **** C **** SET LXX, LYY, LXY C **** DO 14 K=1,KMAX DO 14 I=1,LEN1 V3(I,K,1)=V1(I,K,1)*T1(I)*1.E-1*QE*avo/V6(I,K,1) V4(I,K,1)=V3(I,K,1)*T4(I)**2 V5(I,K,1)=V2(I,K,1)*T1(I)*1.E-1*QE*avo/V6(I,K,1) ! V3(I,K,1) = merge(V3(I,K,1),1.E-20,V3(I,K,1)-1.E-20>=0.) ! V4(I,K,1) = merge(V4(I,K,1),1.E-20,V4(I,K,1)-1.E-20>=0.) ! V5(I,K,1) = merge(V5(I,K,1),1.E-20,V5(I,K,1)-1.E-20>=0.) if (v3(i,k,1) < 1.e-20) v3(i,k,1) = 1.e-20 if (v4(i,k,1) < 1.e-20) v4(i,k,1) = 1.e-20 if (v5(i,k,1) < 1.e-20) v5(i,k,1) = 1.e-20 14 CONTINUE NLXXK=NLXX+1 NLYYK=NLYY+1 NLXYK=NLXY+1 DO 15 I=1,LEN2-LEN1 F(I,NLXXK)=SQRT(V3(I,1,1)*V3(I,2,1)) F(I,NLYYK)=SQRT(V4(I,1,1)*V4(I,2,1)) F(I,NLXYK)=SQRT(V5(I,1,1)*V5(I,2,1)) 15 CONTINUE C **** C **** BOUNDARIES C **** NLXXK=NLXX NLYYK=NLYY NLXYK=NLXY DO 16 I=1,LEN1 F(I,NLXXK)=SQRT(V3(I,1,1)**3/V3(I,2,1)) F(I,NLYYK)=SQRT(V4(I,1,1)**3/V4(I,2,1)) F(I,NLXYK)=SQRT(V5(I,1,1)**3/V5(I,2,1)) F(I,NLXXK+KMAX)=SQRT(V3(I,KMAX,1)**3/V3(I,KMAXM1,1)) F(I,NLYYK+KMAX)=SQRT(V4(I,KMAX,1)**3/V4(I,KMAXM1,1)) F(I,NLXYK+KMAX)=SQRT(V5(I,KMAX,1)**3/V5(I,KMAXM1,1)) 16 CONTINUE NLXYK = NLXY-1 DO 19 K =1,KMAXP1 NLXYK = NLXYK+1 DO 19 I = 1,IMAXP4 F(I,NLXYK) = F(I,NLXYK)*T4(I) 19 CONTINUE C **** C **** COPY SIGMA1, SIGMA2, Z, H, U, V AND W TO COMMON/DYNAMO/ C **** NZK = NJ+NZ-1 NTK = NJ+NT-1 NMSK = NJ+NMS-1 NUK = NJ+NU-1 NVK = NJ+NV-1 NWK = NJ+NW-1 DO 17 K = 1,KMAX NZK = NZK+1 NTK = NTK+1 NMSK = NMSK+1 NUK = NUK+1 NVK = NVK+1 NWK = NWK+1 DO 17 I = 1,IMAX+1 SIGMA1(I,J,K) = V1(I+2,K,1)*1.E10*QE/T1(I+2) SIGMA2(I,J,K) = V2(I+2,K,1)*1.E10*QE/T1(I+2) Z(I,J,K) = F(I+2,NZK) H(I,J,K) = gask*F(I+2,NTK)/(.5*(F(I+2,NMSK)+F(I+2,NMSK+1)) 1 *grav) U(I,J,K) = F(I+2,NUK) V(I,J,K) = F(I+2,NVK) W(I,J,K) = .5*(F(I+2,NWK)+F(I+2,NWK+1))*H(I,J,K) 17 CONTINUE NZK = NZK+1 DO 18 I = 1,IMAX+1 Z(I,J,KMAXP1) = F(I+2,NZK) 18 CONTINUE RETURN END