#include "dims.h" SUBROUTINE LAMDAS use cons_module,only: len1,len2,kmax,kmaxm1,kmaxp1,imax,imaxp4, | rmassinv,expz,gask,grav,avo,p0,boltz,dipmin 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 "trgm.h" real :: bxm,bx,bxp,by,byp,bz,bzp,bmod,bmodp COMMON/MAGFLD/BXM(ZIMXP,2),BX(ZIMXP,ZJMX),BXP(ZIMXP,4), 1 BY(ZIMXP,ZJMX),BYP(ZIMXP,4),BZ(ZIMXP,ZJMX),BZP(ZIMXP,4), 2 BMOD(ZIMXP,ZJMX),BMODP(ZIMXP,2) 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),qe,v4tmp integer :: i,nps1k,nps2k,ntik,ntk,l,m,no2pk,nopk,nnopk,nmsk,k, | ntek,nlxxk,nlyyk,nlxyk,nzk,nuk,nvk,nwk real,dimension(zimxp,zkmxp) :: ! for plotting | ped_plt, hall_plt, z_plt, h_plt, u_plt, v_plt, w_plt ! DATA RMS/32.,16.,30./,QE/1.602E-19/ 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 **** call addfsech('O2P_LAM',' ',' ',f(1,nj+no2p),zimxp,zkmxp,zkmx,j) call addfsech('OP_LAM' ,' ',' ',f(1,nj+nop) ,zimxp,zkmxp,zkmx,j) call addfsech('NOP_LAM',' ',' ',f(1,nj+nnop),zimxp,zkmxp,zkmx,j) DO 1 I=1,LEN1 T1(I)=BMOD(I,J) ! bgauss(lon0:lon1) T2(I)=9.6489E3*T1(I) ! omega_i(lon0:lon1) T3(I)=1.7588028E7*T1(I) ! omega_e(lon0:lon1) ! DIPMIN is in trgm.h, initialized in types.f: if (abs(DIPMAG(I,J)) >= dipmin) then t4(i) = sin(dipmag(i,j)) else t4(i) = sin(dipmin) 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(1) V2(I,1,1)=F(I,NPS2K)*rmassinv(2) V3(I,1,1)=(1.-F(I,NPS1K)-F(I,NPS2K))*rmassinv(3) if (v3(i,1,1) < 1.e-20) v3(i,1,1) = 1.e-20 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 ! call addfsech('VO2',' ',' ',v1,zimxp,zkmxp,zkmx,j) ! call addfsech('VO' ,' ',' ',v2,zimxp,zkmxp,zkmx,j) ! call addfsech('VN2',' ',' ',v3,zimxp,zkmxp,zkmx,j) ! call addfsech('VT' ,' ',' ',v6,zimxp,zkmxp,zkmx,j) C **** C **** RM1 = NEW(L,M) C **** C **** GENERATE NUMERICAL FACTORS CALL NEW(V6,RM1) ! call addfsech('RMO2_1',' ',' ',rm1(1,1,1,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMO2_2',' ',' ',rm1(1,1,2,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMO2_3',' ',' ',rm1(1,1,3,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMO_1' ,' ',' ',rm2(1,1,1,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMO_2' ,' ',' ',rm2(1,1,2,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMO_3' ,' ',' ',rm2(1,1,3,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMN2_1',' ',' ',rm3(1,1,1,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMN2_2',' ',' ',rm3(1,1,2,1),zimxp,zkmxp,zkmx,j) ! call addfsech('RMN2_3',' ',' ',rm3(1,1,3,1),zimxp,zkmxp,zkmx,j) 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 ! call addfsech('VO2',' ',' ',v1,zimxp,zkmxp,zkmx,j) ! call addfsech('VO' ,' ',' ',v2,zimxp,zkmxp,zkmx,j) ! call addfsech('VN2',' ',' ',v3,zimxp,zkmxp,zkmx,j) 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 ! call addfsech('XNMBARM',' ',' ',v6,zimxp,zkmxp,zkmx,j) ! call addfsech('VO2',' ',' ',v1,zimxp,zkmxp,zkmx,j) ! call addfsech('VO' ,' ',' ',v2,zimxp,zkmxp,zkmx,j) ! call addfsech('VN2',' ',' ',v3,zimxp,zkmxp,zkmx,j) 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)) if ((1.-F(I,NPS1K)-F(I,NPS2K))*rmassinv(3) >= 1.E-20) then v4tmp =(1.-F(I,NPS1K)-F(I,NPS2K))*rmassinv(3) else v4tmp = 1.e-20 endif V4(I,1,1)=V6(I,1,1)*(2.33E-11*v4tmp* | F(I,NTEK)*(1.-1.21E-4*F(I,NTEK))+1.82E-10*F(I,NPS1K)* | rmassinv(1)*V5(I,1,1)*(1.+3.6E-2*V5(I,1,1))+8.9E-11* | F(I,NPS2K)*rmassinv(2)*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 ! call addfsech('NU_OMEGA',' ',' ',v4,zimxp,zkmxp,zkmx,j) 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 ! call addfsech('PED_1',' ',' ',rm2(1,1,1,1),zimxp,zkmxp,zkmx,j) ! call addfsech('PED_2',' ',' ',rm2(1,1,2,1),zimxp,zkmxp,zkmx,j) ! call addfsech('PED_3',' ',' ',rm2(1,1,3,1),zimxp,zkmxp,zkmx,j) ! call addfsech('HALL_1',' ',' ',rm3(1,1,1,1),zimxp,zkmxp,zkmx,j) ! call addfsech('HALL_2',' ',' ',rm3(1,1,2,1),zimxp,zkmxp,zkmx,j) ! call addfsech('HALL_3',' ',' ',rm3(1,1,3,1),zimxp,zkmxp,zkmx,j) 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 ! call addfsech('XNMBARM',' ',' ',v6,zimxp,zkmxp,zkmx,j) ! call addfsech('PEDXIONS',' ',' ',v1,zimxp,zkmxp,zkmx,j) ! call addfsech('HALXIONS',' ',' ',v2,zimxp,zkmxp,zkmx,j) 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)*T5(I) V5(I,K,1)=V2(I,K,1)*T1(I)*1.E-1*QE*avo/V6(I,K,1) 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)) C **** C **** C **** LEVELS 1 AND KMAXP1 C F(I,NLXXK)= 1.5*V3(I,1,1)-0.5*V3(I,2,1) C F(I,NLYYK)= 1.5*V4(I,1,1)-0.5*V4(I,2,1) C F(I,NLXYK)= 1.5*V5(I,1,1)-0.5*V5(I,2,1) C F(I,NLXXK+KMAX)=1.5*V3(I,KMAX,1)-0.5*V3(I,KMAXM1,1) C F(I,NLYYK+KMAX)=1.5*V4(I,KMAX,1)-0.5*V4(I,KMAXM1,1) C F(I,NLXYK+KMAX)=1.5*V5(I,KMAX,1)-0.5*V5(I,KMAXM1,1) C **** C **** 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 ! call addfsech('LXX',' ',' ',f(1,nlxx),zimxp,zkmxp,zkmx,j) ! call addfsech('LYY',' ',' ',f(1,nlyy),zimxp,zkmxp,zkmx,j) ! call addfsech('LXY',' ',' ',f(1,nlxy),zimxp,zkmxp,zkmx,j) ! call addfsech('LYX',' ',' ',f(1,nlyx),zimxp,zkmxp,zkmx,j) 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) ped_plt(i+2,k) = sigma1(i,j,k) hall_plt(i+2,k) = sigma2(i,j,k) z_plt(i+2,k) = z(i,j,k) h_plt(i+2,k) = h(i,j,k) u_plt(i+2,k) = u(i,j,k) v_plt(i+2,k) = v(i,j,k) w_plt(i+2,k) = w(i,j,k) 17 CONTINUE NZK = NZK+1 DO 18 I = 1,IMAX+1 Z(I,J,KMAXP1) = F(I+2,NZK) z_plt(i+2,kmaxp1) = z(i,j,kmaxp1) 18 CONTINUE ! These calls are in dynamo.F in tiegcm1: ! call addfsech('PED' ,' ',' ',ped_plt ,zimxp,zkmxp,zkmx,j) ! call addfsech('HALL' ,' ',' ',hall_plt,zimxp,zkmxp,zkmx,j) ! call addfsech('SIGMAPED',' ',' ',ped_plt ,zimxp,zkmxp,zkmx,j) ! call addfsech('SIGMAHAL',' ',' ',hall_plt,zimxp,zkmxp,zkmx,j) ! call addfsech('ZPOTEN' ,' ',' ',z_plt ,zimxp,zkmxp,zkmx,j) ! call addfsech('SCHEIGHT',' ',' ',h_plt ,zimxp,zkmxp,zkmx,j) ! call addfsech('UNVEL' ,' ',' ',u_plt ,zimxp,zkmxp,zkmx,j) ! call addfsech('VNVEL' ,' ',' ',v_plt ,zimxp,zkmxp,zkmx,j) ! call addfsech('WNVEL' ,' ',' ',w_plt ,zimxp,zkmxp,zkmx,j) RETURN END C