SUBROUTINE QRJ CDIR$ VFUNCTION EXPHF implicit none C **** CALCULATE HEATING AND DISSOCIATION RATES C **** CHAPMAN INTEGRAL FOR EACH SPECIES AT EACH GRID C **** POINT IN NNO2, NNO, NNN2 (CALCULATED BY SUB CHAPMN) include "params.h" include "blnk.h" include "vscr.h" include "cons.h" include "index.h" include "buff.h" include "phys.h" include "crates.h" include "cmpbnd.h" real :: rmn4s,rmn2d,rmno,brn2d,cee COMMON/MASS/ RMN4S,RMN2D,RMNO,BRN2D,CEE integer,parameter :: lmax=59 real :: EUVEFF(ZKMXP),SIGEUV(3,lmax),FEUV(lmax),RLMEUV(lmax), | FSRC(15),SIGSRC(15),RLMSRC(15),SIGMAS(6,lmax),QUENCH(4), | SFLUX(lmax),BRN2(lmax),BRO2(lmax) common/qrj_coeff/ euveff,sigeuv,feuv,rlmeuv,fsrc,sigsrc,rlmsrc, | sigmas,quench,sflux,brn2,bro2 ! ! Local: integer :: n,m,nno2k,nnok,nnn2k,i,nps1k,nps2k, | npn4sk,nqk,nrjk,nqtefk,nqop2pk,nqop2dk,k,nnk,ntk,nmsk,nqo2pk, | nqopk,nqn2pk,nqnpk,nqnopk,npnok real :: do2,ad,bd,aband,do22,e3,f107,f107a,factor,bband,cband real :: exphf ! DO2=8.203E-12 AD=5.435E-9 BD=1.232E-6 ABAND=0.143 BBAND=9.64E8 CBAND=9.03E-19 DO22=1.1407E-11 E3=0.33 C **** C **** C **** COMPUTE S14 = RSQ, S15 = RSP C **** C **** C **** S1 = TAU(R), S2 = TAU(Q) C **** NNO2K = NNO2 NNOK = NNO NNN2K = NNN2 DO 6 I = 1,LEN3 S1(I,1) = SIGEUV(1,49)*F(I,NNO2K)+SIGEUV(2,49)*F(I,NNOK)+ 1 SIGEUV(3,49)*F(I,NNN2K) S2(I,1) = SIGEUV(1,20)*F(I,NNO2K)+SIGEUV(2,20)*F(I,NNOK)+ 1 SIGEUV(3,20)*F(I,NNN2K) C **** C **** S3 = TAU(1), S4 = TAU(2), S5 = TAU(3) C **** S3(I,1) = 1.3*S1(I,1) S4(I,1) = 2.0*S1(I,1) S5(I,1) = 2.5*S1(I,1) C **** C **** IF(TAU.GT.9.0) TAU = 9.0 C **** S1(I,1) = merge(S1(I,1),9.,9.-S1(I,1)>=0.) S2(I,1) = merge(S2(I,1),9.,9.-S2(I,1)>=0.) S3(I,1) = merge(S3(I,1),9.,9.-S3(I,1)>=0.) S4(I,1) = merge(S4(I,1),9.,9.-S4(I,1)>=0.) S5(I,1) = merge(S5(I,1),9.,9.-S5(I,1)>=0.) C **** C **** TAU(N) = EXP(-TAU(N)) FOR N = 1,3,1 C **** S3(I,1) = exphf(-S3(I,1)) S4(I,1) = exphf(-S4(I,1)) S5(I,1) = exphf(-S5(I,1)) C **** C **** S6 = EXP(-TAU(R)), S7 = EXP(-TAU(Q)) C **** S6(I,1) = exphf(-S1(I,1)) S7(I,1) = exphf(-S2(I,1)) C **** C **** S14 = RSP, S15 = RSQ C **** S14(I,1) = S6(I,1)+2.*(S3(I,1)+S4(I,1)+S5(I,1)) S15(I,1) = 1.5*S6(I,1)/(S14(I,1)+S2(I,1)/S1(I,1)*S7(I,1)) S14(I,1) = 2.4*S6(I,1)/S14(I,1) 6 CONTINUE C **** C **** C **** S1 = PSI(O2), S2 = PSI(O), S3 = PSI(N2), S4 = PSI(N4S) C **** C **** NPS1K=NJ+NPS NPS2K=NJ+NPS2 NPN4SK=NJ+NPN4S C **** C **** LEVELS 2 THRU KMAXP1 C **** DO 7 I=1,LEN2 S1(I,2)=.5*(F(I,NPS1K)+F(I,NPS1K+1)) S2(I,2)=.5*(F(I,NPS2K)+F(I,NPS2K+1)) C S4(I,2)=.5*(F(I,NPN4SK)+F(I,NPN4SK+1)) S4(I,2)=0. 7 CONTINUE C **** C **** LEVEL 1 C **** DO 8 I=1,LEN1 S1(I,1) = .5*((B(I,1,1)+1.)*F(I,NPS1K)+B(I,1,2)*F(I,NPS2K)+ 1 FB(I,1)) S2(I,1) = .5*(B(I,2,1)*F(I,NPS1K)+(B(I,2,2)+1.)*F(I,NPS2K)+ 1 FB(I,2)) C S4(I,1) = F(I,NPN4SK)**1.5/SQRT(F(I,NPN4SK+1)) S4(I,1) = 0. 8 CONTINUE C **** C **** NOW CALCULATE S3 = PSI(N2) C **** DO 9 I=1,LEN3 S3(I,1)=1.-S1(I,1)-S2(I,1) 9 CONTINUE C **** C **** C **** CALCULATE IONIZATION CONTRIBUTIONS IN S8 - S13 C **** C **** C **** S8 = O2 IONIZATION C **** S9 = O IONIZATION C **** S10 = N2 IONIZATION C **** S11 = N4S IONIZATION C **** S12 = O2 IONIZATION * BRO2 C **** S13 = N2 IONIZATION * BRN2 C **** C **** C **** CONTRIBUTIONS TO Q, J AND QTEF FROM EUV C **** C **** INITIALIZE Q, RJ AND QTEF C **** NQK = NQ NRJK = NRJ NQTEFK = NQTEF NQOP2PK = NQOP2P NQOP2DK = NQOP2D DO 10 I = 1,LEN3 F(I,NQK) = 0. F(I,NRJK) = 0. F(I,NQTEFK) = 0. F(I,NQOP2PK) = 0. F(I,NQOP2DK) = 0. 10 CONTINUE C **** C **** INITIALIZE S8 - S13 C **** DO 11 I=1,6*LEN3 S13(I,1)=0. 11 CONTINUE C **** C **** SUMMATION OVER WAVE LENGTH C **** DO 12 N=16,lmax C **** C **** S5=SUM(O2,O,N2)(SIGMA*CHAPMAN) C **** S6=SUM(O2,O,N2)(SIGMA*PSI/RMASS) C **** S7=SUM(O2,O,N2,N4S)(SIGMAS) C **** C **** INITIALIZE S5, S6, S7 C **** DO 13 I=1,3*LEN3 S7(I,1)=0. 13 CONTINUE C **** C **** SUMMATION OVER O2, O, N2 C **** DO 14 M=1,3 K=(3-M)*KMAXP1+1 NNK=NNO2+(M-1)*KMAXP1 DO 14 I=1,LEN3 S5(I,1)=S5(I,1)+SIGEUV(M,N)*F(I,NNK) S6(I,1)=S6(I,1)+SIGEUV(M,N)*S3(I,K)/RMASS(M) S7(I,1)=S7(I,1)+SIGMAS(M,N) 14 CONTINUE DO 15 I = 1,LEN3 C **** C **** ADD CONTRIBUTIONS FROM N4S TO S7 C **** S7(I,1) = S7(I,1)+SIGMAS(4,N) C **** C **** S5 = F(LAMDA)*EXP(-TAU) C **** S5(I,1) = FEUV(N)*exphf(-S5(I,1)) C **** C **** ADD CONTRIBUTIONS TO Q AND RJE C **** F(I,NQK) = F(I,NQK)+C(60)/RLMEUV(N)*S5(I,1)*S6(I,1) F(I,NRJK) = F(I,NRJK)+S5(I,1)*S7(I,1) F(I,NQOP2PK) = F(I,NQOP2PK)+S5(I,1)*SIGMAS(6,N)* 1 S2(I,1)/RMASS(2) F(I,NQOP2DK) = F(I,NQOP2DK)+S5(I,1)*SIGMAS(5,N)* 1 S2(I,1)/RMASS(2) 15 CONTINUE C **** C **** ADD IONIZATION FOR O2, O, N2 TO S8, S9, S10 C **** DO 16 M=1,3 K=(3-M)*KMAXP1+1 DO 16 I=1,LEN3 S10(I,K)=S10(I,K)+SIGMAS(M,N)*S5(I,1)*S3(I,K)/RMASS(M) 16 CONTINUE DO 17 I = 1,LEN3 C **** C **** ADD EFFECTIVE N2 IONIZATION TO QTEF C **** F(I,NQTEFK) = F(I,NQTEFK)+(SIGEUV(3,N)-SIGMAS(3,N))*S5(I,1)* 1 S3(I,1)/RMASS(3)*2. C **** C **** ADD N4S IONIZATION TO S11 C **** S11(I,1) = S11(I,1)+SIGMAS(4,N)*S5(I,1)*S4(I,1)/RMN4S C **** C **** ADD IONIZATION CONTRIBUTIONS TO S12 AND S13 C **** S12(I,1) = S12(I,1)+SIGMAS(1,N)*S5(I,1)*S1(I,1)/RMASS(1)* 1 BRO2(N) S13(I,1) = S13(I,1)+SIGMAS(3,N)*S5(I,1)*S3(I,1)/RMASS(3)* 1 BRN2(N) 17 CONTINUE 12 CONTINUE C **** C **** MULTIPLY Q BY EFICIENCY FACTOR C **** C **** NQK=NQK-1 DO 18 K=1,KMAXP1 NQK=NQK+1 DO 18 I=1,LEN1 F(I,NQK)=F(I,NQK)*EUVEFF(K)*C(85) 18 CONTINUE C **** C **** C **** CALCULATE N*MBAR IN S5 C **** C **** C **** S6 = T(TOTAL,K) C **** C **** LEVEL 1 C **** NTK=NJ+NT+KMAX DO 19 I=1,LEN1 S6(I,1)=F(I,NTK)+T0(1) 19 CONTINUE C **** C **** LEVELS 2 THRU KMAX C **** NTK=NJ+NT-1 DO 20 K=2,KMAX NTK=NTK+1 DO 20 I=1,LEN1 S6(I,K)=.5*(F(I,NTK)+F(I,NTK+1))+T0(K) 20 CONTINUE C **** C **** LEVEL KMAXP1 C **** NTK=NJ+NT+KMAX-1 DO 21 I=1,LEN1 S6(I,KMAXP1)=F(I,NTK)+T0(KMAXP1) 21 CONTINUE C **** C **** S7 = EXPS (K) C **** DO 22 K=1,KMAX DO 22 I=1,LEN1 S7(I,K)=C(87)*EXPS(K) 22 CONTINUE DO 23 I=1,LEN1 S7(I,KMAXP1)=C(86)*EXPS(KMAX) 23 CONTINUE C **** C **** S5 = N*MBAR (K) C **** NMSK=NJ+NMS DO 24 I=1,LEN3 S5(I,1)=C(81)*S7(I,1)*F(I,NMSK)/(C(84)*S6(I,1)) 24 CONTINUE C **** C **** C **** CALCULATE CONTRIBUTIONS TO NQO2P, NQOP, NQN2P, NQNP C **** AND NQTEF C **** C **** NQO2PK = NQO2P NQOPK = NQOP NQN2PK = NQN2P NQNPK = NQNP NQTEFK = NQTEF NQOP2PK = NQOP2P NQOP2DK = NQOP2D DO 25 I=1,LEN3 F(I,NQO2PK) = F(I,NQO2PK)+(S8(I,1)*(1.+S15(I,1))-S12(I,1))* 1 S5(I,1) F(I,NQN2PK) = F(I,NQN2PK)+(S10(I,1)*(1.+S14(I,1))-S13(I,1))* 1 S5(I,1) DISN2P(I,1) = S10(I,1)*S14(I,1)*S5(I,1) F(I,NQNPK) = F(I,NQNPK)+(S11(I,1)+S13(I,1))*S5(I,1) F(I,NQTEFK) = F(I,NQTEFK)*S5(I,1) S4(I,1) = F(I,NQOP2PK)+F(I,NQOP2DK)+S9(I,1)+S12(I,1) F(I,NQOP2PK) = (F(I,NQOP2PK)+S4(I,1)*S14(I,1)*0.22)*S5(I,1) F(I,NQOP2DK) = (F(I,NQOP2DK)+S4(I,1)*S14(I,1)*0.24)*S5(I,1) F(I,NQOPK) = F(I,NQOPK)+(S9(I,1)+S12(I,1)+S4(I,1)*S14(I,1) 1 *0.56)*S5(I,1) 25 CONTINUE C **** C **** CALCULATE NO IONIZATION AND ADD TO NQNOP C **** NQNOPK=NQNOP NNO2K=NNO2 NPNOK=NJ+NPNO C **** LEVEL 1 DO 26 I=1,LEN1 F(I,NQNOPK)=F(I,NQNOPK)+BETA9(I,1)*F(I,NPNOK) 1 *S5(I,1)/RMNO 26 CONTINUE C **** LEVELS 2 THRU KMAXP1 DO 27 I=LEN1+1,LEN3 F(I,NQNOPK)=F(I,NQNOPK)+BETA9(I,1)*.5*(F(I,NPNOK)+F(I,NPNOK-1)) 1 *S5(I,1)/RMNO 27 CONTINUE DO 28 K=1,KMAXP1 FACTOR=C(85)*C(81)/C(57)*EXPS(1)*C(86)**(2*K-3) DO 28 I=1,LEN1 S8(I,K)=FACTOR/((S1(I,K)/RMASS(1)+S2(I,K)/RMASS(2)+ 1 S3(I,K)/RMASS(3))*S6(I,K)) 28 CONTINUE DO 29 I=1,LEN3 S8(I,1)=S8(I,1)*(QUENCH(1)*S3(I,1)/RMASS(3)+ 1 QUENCH(2)*S1(I,1)/RMASS(1)) S8(I,1)=QUENCH(3)*S8(I,1)/(QUENCH(4)+S8(I,1)) 29 CONTINUE NQK=NQ NRJK=NRJ C **** C **** S7=SUM OVER WAVE LENGTH(SIGMA*F*EXP(-SIGMA*CHAPMAN)* C **** (HC/LAMDA-DO2)) C **** INITIALIZE S7 C **** DO 30 I=1,LEN3 S7(I,1)=0. 30 CONTINUE C **** C **** SUMMATION OVER WAVE LENGTH C **** ! ! If there is a crash here, check that init_sigmas has been called ! by start, and init_sflux has been called by advnce. ! NNO2K=NNO2 DO 31 N=1,15 DO 31 I=1,LEN3 S9(I,1)=SIGSRC(N)*FSRC(N)*exphf(-SIGSRC(N)*F(I,NNO2K)) S7(I,1)=S7(I,1)+S9(I,1)* 1 (C(60)/RLMSRC(N)-DO22+S8(I,1)) C **** C **** UPDATE RJ C **** F(I,NRJK)=F(I,NRJK)+S9(I,1) 31 CONTINUE C **** C **** UPDATE Q C **** DO 32 I=1,LEN3 F(I,NQK)=F(I,NQK)+S7(I,1)*C(85)*S1(I,1)/RMASS(1) 32 CONTINUE C **** CONTRIBUTIONS FROM SCHUMANN RUNGE BANDS C **** S7=P3*F NNO2K=NNO2 DO 33 I=1,LEN3 S8(I,1)=SQRT(F(I,NNO2K)) S7(i,1) = merge(1./(ABAND*F(I,NNO2K)+BBAND*S8(I,1)),CBAND, + F(I,NNO2K)-1.E18>=0.0)*(1.+0.11*(C(61)-65.)/165.) F(I,NQK)=F(I,NQK)+S7(I,1)*C(85)*S1(I,1)/RMASS(1)*E3 F(I,NRJK)=F(I,NRJK)+S7(I,1)/DO2 33 CONTINUE RETURN END !--------------------------------------------------------------- subroutine init_sflux implicit none ! ! Flux initialization once per time step, called from advnce. ! (this call moved from qrj in multi-task version) ! include "params.h" include "cons.h" integer,parameter :: lmax=59 ! ! Local: real :: EUVEFF(ZKMXP),SIGEUV(3,lmax),FEUV(lmax),RLMEUV(lmax), | FSRC(15),SIGSRC(15),RLMSRC(15),SIGMAS(6,lmax),QUENCH(4), | SFLUX(lmax),BRN2(lmax),BRO2(lmax) real :: EUVFLX(37),wave1(lmax),wave2(lmax) ! output from euvac common/qrj_coeff/ euveff,sigeuv,feuv,rlmeuv,fsrc,sigsrc,rlmsrc, | sigmas,quench,sflux,brn2,bro2 real :: flya,f107,f107a,hlybr,fexvir,hlya,heiew,xuvfac integer :: iscale,n,nn F107 = C(61) F107A = C(62) ISCALE = 0 HLYBR = 0. FEXVIR = 0. HLYA = 3.E+11+0.4E+10*(C(61)-70.) HEIEW = 0. C XUVFAC =0. C XUVFAC = 2.0 - (C(61)-68.0) / (243.0-68.0) XUVFAC = 4.0 - (C(61)-68.0) / (243.0-68.0) IF (XUVFAC .LT. 1.0) XUVFAC = 1.0 C C ISCALE =0 for Hinteregger contrast ratio method C =1 for Hinteregger linear interpolation C =2 for Tobiska EUV91 model C =3 for Woods & Rottman 10 Nov. 1988 measurement C =4 for Woods & Rottman 20 Jun. 1989 measurement C F107 daily 10.7 cm flux (1.E-22 W m-2 Hz-1) C F107A 81-day centered average 10.7 cm flux C HLYBR ratio of H Ly-b 1026A flux to solar minimum value (optional) C FEXVIR ratio of Fe XVI 335A flux to solar minimum value (optional) C HLYA H Lyman-alpha flux (photons cm-2 s-1) (optional) C HEIEW He I 10830A equivalent width, (milliAngstroms) (optional) C XUVFAC factor for scaling flux 16-250A (optional) C WAVE1 longwave bound of spectral intervals (Angstroms) C WAVE2 shortwave bound of intervals (= WAVE1 for indiv. lines) C SFLUX scaled solar flux returned by subroutine (photons cm-2 s-1) C CALL SSFLUX(ISCALE,F107,F107A,HLYBR,FEXVIR,HLYA, | HEIEW, XUVFAC, WAVE1, WAVE2, SFLUX) CALL EUVAC(F107,F107A,EUVFLX) C **** C **** TRANSFER VALUES OF SFLUX TO APPROPRIATE SLOTS C **** DO 4 N = 16,lmax FEUV(N) = SFLUX(N) 4 CONTINUE DO 5 N = 1,15 FSRC(N) = SFLUX(N) 5 CONTINUE DO 55 N=1,37 NN = N+15 FEUV(NN) = EUVFLX(N) 55 CONTINUE return end !--------------------------------------------------------------- block data set_qrj_coeff include "params.h" integer,parameter :: lmax=59 real :: EUVEFF(ZKMXP),SIGEUV(3,lmax),FEUV(lmax),RLMEUV(lmax), | FSRC(15),SIGSRC(15),RLMSRC(15),SIGMAS(6,lmax),QUENCH(4), | SFLUX(lmax),BRN2(lmax),BRO2(lmax) common/qrj_coeff/ euveff,sigeuv,feuv,rlmeuv,fsrc,sigsrc,rlmsrc, | sigmas,quench,sflux,brn2,bro2 C **** C **** TOTAL ABSORPTION COEFFICIENTS C **** C **** SIGEUV(1,N) = SIGA(O2) C **** SIGEUV(2,N) = SIGA(O) C **** SIGEUV(3,N) = SIGA(N2) C **** C **** TOTAL IONIZATION COEFFICIENTS C **** C **** SIGMAS(1,N) = SIGI(O2) C **** SIGMAS(2,N) = SIGI(O+(4S)) C **** SIGMAS(3,N) = SIGI(N2) C **** SIGMAS(4,N) = SIGI(N) C **** SIGMAS(5,N) = SIGI(O+(2D)) C **** SIGMAS(6,N) = SIGI(O+(2P)) C **** C **** C DATA RLMEUV/0.17250E-04, 0.16750E-04, 0.16250E-04, 0.15750E-04, + 0.15250E-04, 0.14750E-04, 0.14250E-04, 0.13750E-04, + 0.13250E-04, 0.12750E-04, 0.12250E-04, 0.12157E-04, + 0.11750E-04, 0.11250E-04, 0.10750E-04, 0.10250E-04, + 0.10319E-04, 0.10257E-04, 0.97500E-05, 0.97702E-05, + 0.92500E-05, 0.87500E-05, 0.82500E-05, 0.77500E-05, + 0.78936E-05, 0.77041E-05, 0.76515E-05, 0.72500E-05, + 0.70331E-05, 0.67500E-05, 0.62500E-05, 0.62973E-05, + 0.60976E-05, 0.57500E-05, 0.58433E-05, 0.55437E-05, + 0.52500E-05, 0.47500E-05, 0.46522E-05, 0.42500E-05, + 0.37500E-05, 0.36807E-05, 0.32500E-05, 0.30378E-05, + 0.30331E-05, 0.27500E-05, 0.28415E-05, 0.25630E-05, + 0.22500E-05, 0.17500E-05, 0.12500E-05, 0.75000E-06, + 0.41000E-06, 0.27500E-06, 0.19500E-06, 0.12000E-06, + 0.60000E-07, 0.30000E-07, 0.15000E-07/ DATA(SIGEUV(1,N),N=1,lmax)/ + 0.50, 1.50, 3.40, 6.00,10.00,13.00, + 15.00,12.00, 2.20, 0.30, 3.00, 0.01, + 0.30, 0.10, 1.00, 1.10, 1.00, 1.60, + 16.53, 4.00,15.54, 9.85,20.87,27.09, + 26.66,25.18,21.96,29.05,25.00,26.27, + 26.02,25.80,26.10,25.04,22.00,25.59, + 24.06,21.59,20.40,19.39,18.17,18.40, + 17.19,16.80,16.80,15.10,15.70,13.20, + 10.60, 7.10, 4.00, 1.18, 0.32, 0.10, + 1.02, 0.14, .024, .004, .0004/ DATA(SIGEUV(2,N),N=1,lmax)/ + 18 * 0.00, + 0.00, 0.00, 2.12, 4.18, 4.38, 4.23, + 4.28, 4.18, 4.18, 8.00,11.35,10.04, + 12.21,12.22,12.23,11.90,12.17,12.13, + 11.91,11.64,11.25,11.21, 9.64, 9.95, + 8.67, 7.70, 7.68, 6.61, 7.13, 6.05, + 5.30, 2.90, 1.60, 0.59, 0.16, 0.05, + 0.51, 0.07, .012, .002, .0002/ DATA(SIGEUV(3,N),N=1,lmax)/ + 18 * 0.00, + 36.16, 0.70,16.99,46.63,15.05,30.71, + 19.26,26.88,35.46,30.94,26.30,29.75, + 23.22,23.20,23.10,22.38,23.20,24.69, + 24.53,21.85,21.80,21.07,17.51,18.00, + 13.00,11.60,11.60,10.30,10.60, 9.70, + 8.00, 4.40, 1.90, 0.60, 0.24, 1.16, + 0.48, 0.09, .015, .003, .0003/ DATA(SIGMAS(1,N),N=1,lmax)/ + 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, + 0.00, 0.00, 0.00, 0.00, 0.00, 0.00, + 0.00, 0.00, 0.00, 0.27, 0.00, 1.00, + 12.22, 2.50, 9.34, 4.69, 6.12, 9.39, + 11.05, 9.69, 8.59,23.81,23.00,22.05, + 25.94,25.80,26.10,25.04,22.00,25.59, + 24.06,21.59,20.40,19.39,18.17,18.40, + 17.19,16.80,16.80,15.10,15.70,13.20, + 10.60, 7.10, 4.00, 1.18, 0.32, 0.10, + 1.02, 0.14, .024, .004, .0004/ DATA(SIGMAS(2,N),N=1,lmax)/ + 20 * 0.00, 2.12, 4.18, 4.38, 4.23, + 4.28, 4.18, 4.18, 4.20, 4.91, 4.01, + 3.78, 3.79, 3.67, 3.45, 3.53, 3.52, + 3.45, 3.26, 3.15, 3.03, 2.51, 2.59, + 2.25, 1.93, 1.92, 1.65, 1.78, 1.51, + 1.38, 0.78, 0.46, 0.18, 0.05, 0.015, + 0.015, 0.02, 0.004, 0.0006, 0.00006/ DATA(SIGMAS(3,N),N=1,lmax)/ + 18 * 0.00, + 0.00, 0.00, 0.00, 0.00, 0.00,16.75, + 10.18,18.39,23.77,23.20,23.00,25.06, + 23.22,23.20,23.10,22.38,23.20,24.69, + 24.53,21.85,21.80,21.07,17.51,18.00, + 13.00,11.60,11.60,10.30,10.60, 9.70, + 8.00, 4.40, 1.90, 0.60, 0.24, 1.16, + 0.48, 0.09, .015, .003, .0003/ DATA(SIGMAS(4,N),N=1,lmax)/ + 22 * 0.00, + 1.00, 1.00, 1.00, 1.00, 1.10, 1.10, + 1.10, 1.20, 1.20, 1.20, 1.20, 1.20, + 1.10, 1.20, 1.15, 1.10, 1.00, 1.00, + 1.00, 0.70, 0.80, 0.65, 0.60, 0.60, + 0.50, 0.50, 0.40, 0.35, 0.25, 0.20, + 0.10, 0.10, 0.10, 0.05, 0.01, + 2*0.0/ DATA(SIGMAS(5,N),N=1,lmax)/ + 27 * 0.00, 3.80, 6.44, 5.52, 5.49, + 5.50, 5.50, 5.36, 5.48, 5.46, 5.36, + 5.24, 5.06, 4.71, 3.86, 3.98, 3.47, + 2.85, 2.84, 2.38, 2.64, 2.12, 1.86, + 0.99, 0.51, 0.19, 0.05, 0.015, 0.015, + 0.02, 0.004, 0.0006, 0.00006/ DATA(SIGMAS(6,N),N=1,lmax)/ + 29 * 0.00, 0.50, 2.93, 2.93, 3.06, + 3.09, 3.16, 3.15, 3.10, 3.14, 3.04, + 3.48, 3.28, 3.38, 2.95, 2.93, 2.92, + 2.58, 2.71, 2.42, 2.07, 1.13, 0.62, + 0.22, 0.06, 0.02, 0.02, 0.03, 0.004, + 0.0007, 0.00007/ DATA(BRN2(N),N=1,lmax)/ + 36 * 0.00, + 0.01, 0.04, 0.04, 0.03, 0.05, 0.05, + 0.15, 0.20, 0.20, 0.25, 0.32, 0.34, + 11 * 0.36/ DATA(BRO2(N),N=1,lmax)/ + 30 * 0.00, + .025, .036, .045, .120, .155, .189, + .230, 0.20, 0.20, 0.20, 0.23, 0.25, + 0.29, 16 * 0.33/ DATA EUVEFF/29*0.05/ C DATA EUVEFF/29*0.10/ DATA QUENCH/7.E-11,5.E-11,3.1401E-12,9.1E-3/ end ! block data set_qrj_coeff !--------------------------------------------------------------- subroutine init_sigmas implicit none C **** C **** INITIALIZATION UPON FIRST ENTRY C **** ! Called once per run from start. ! (this call was moved from qrj in multi-task version) ! SIGxxx were init by block data set_qrj_coeff (qrj.f) ! include "params.h" integer,parameter :: lmax=59 real :: EUVEFF(ZKMXP),SIGEUV(3,lmax),FEUV(lmax),RLMEUV(lmax), | FSRC(15),SIGSRC(15),RLMSRC(15),SIGMAS(6,lmax),QUENCH(4), | SFLUX(lmax),BRN2(lmax),BRO2(lmax) common/qrj_coeff/ euveff,sigeuv,feuv,rlmeuv,fsrc,sigsrc,rlmsrc, | sigmas,quench,sflux,brn2,bro2 real :: wleuv1,wleuv2,sigao,sigao2,sigan2,sigio,sigio2, | sigin2,brop4s,brop2d,brop2p,sigop2p,sigop2d,sigop4s, | sigin,brn2np,bro2op COMMON/EUV/ WLEUV1(37),WLEUV2(37),SIGAO(37),SIGAO2(37), 1 SIGAN2(37),SIGIO(37),SIGIO2(37),SIGIN2(37), 2 BROP4S(37),BROP2D(37),BROP2P(37), 3 SIGOP2P(37),SIGOP2D(37),SIGOP4S(37), 4 SIGIN(37),BRN2NP(37),BRO2OP(37) ! ! Local: integer :: m,n,nn ! DO 1 M=1,3 DO 1 N=1,lmax SIGEUV(M,N)=SIGEUV(M,N)*1.E-18 SIGMAS(M,N)=SIGMAS(M,N)*1.E-18 1 CONTINUE DO 2 N = 1,lmax SIGMAS(4,N) = SIGMAS(4,N)*1.E-17 SIGMAS(5,N) = SIGMAS(5,N)*1.E-18 SIGMAS(6,N) = SIGMAS(6,N)*1.E-18 2 CONTINUE DO 3 N = 1,15 RLMSRC(N) = RLMEUV(N) SIGSRC(N) = SIGEUV(1,N) 3 CONTINUE ! ! DO 51 loop moved here from original euvac. DO 51 N=1,37 BROP2P(N) = 0. IF(N.GT.14) BROP2P(N) = 1.-BROP2D(N)-BROP4S(N) SIGOP2P(N)=SIGIO(N)*BROP2P(N) SIGOP2D(N)=SIGIO(N)*BROP2D(N) SIGOP4S(N)=SIGIO(N)*BROP4S(N) 51 CONTINUE DO 56 N = 1,37 NN = N+15 ! 16:52 SIGEUV(1,NN) = SIGAO2(N) SIGEUV(2,NN) = SIGAO(N) SIGEUV(3,NN) = SIGAN2(N) SIGMAS(1,NN) = SIGIO2(N) SIGMAS(2,NN) = SIGOP4S(N) SIGMAS(3,NN) = SIGIN2(N) SIGMAS(4,NN) = SIGIN(N) SIGMAS(5,NN) = SIGOP2D(N) SIGMAS(6,NN) = SIGOP2P(N) BRN2(NN) = BRN2NP(N) BRO2(NN) = BRO2OP(N) 56 CONTINUE return end