! SUBROUTINE QRJ C 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) C use input_module,only: f107 use init_module,only: sfeps implicit none include "params.h" include "fcom.h" include "cons.h" include "index.h" include "buff.h" include "phys.h" include "crates.h" include "cmpbnd.h" include "mwt.h" include "qrj_coeff.h" ! ! Local: integer :: n,m,nno2k,nnok,nnn2k,i,nps1k,nps2k,nqk,nrjk,nqtefk, | nqop2pk,nqop2dk,k,nnk,ntk,nmsk,nqo2pk,nqopk,nqn2pk,nqnpk, | nqnopk,npnok,kk real :: rlmeuvinv(lmax),rlmsrcinv(l1-1),rmn4sinv,rmnoinv,factor real,dimension(zkmxp,zimxp) :: | s1k ,s2k ,s3k ,s4k ,s5k ,s6k ,s7k ,s8k ,s9k ,s10k , | s11k ,s12k ,s13k ,s14k ,s15k real,dimension(zkmxp,zimxp) :: | fnno2 ,fnno ,fnnn2 ,fnq ,fnrj ,fnqop2p ,fnqop2d, | fnqtef ,fnt ,fnqo2p ,fnqop ,fnqnp ,fnqn2p ,fnms , | fnqnop ,fnpno ,fnps1 ,fnps2 ,beta9ki ,disn2pki real,dimension(l1:lmax,zkmxp,zimxp) :: s5nk, s6nk, s7nk real,parameter :: | do2=8.203E-12 , | do22=1.1407E-11 , | ad=5.435E-9 , ! not used | bd=1.232E-6 , ! not used | aband=0.143 , ! shumann-runge | bband=9.64E8 , ! shumann-runge | cband=9.03E-19 , ! shumann-runge | e3=0.33 , | hc = 1.9845E-16 ! C(60) ! ! Exec: ! write(6,"('qrj: j=',i2,' number of daytime indices = ',i2, ! | ' number of night indices = ',i2)") j,count(idn==1), ! | count(idn==0) ! For experiment, set idn(:) all daytime: ! idn(:) = 1 ! ! More efficient to multiply by inverse than to divide: ! (see also rmassinv) do i=1,lmax rlmeuvinv(i) = 1./rlmeuv(i) enddo do i=1,l1-1 ! 1,15 rlmsrcinv(i) = 1./rlmsrc(i) enddo rmn4sinv = 1./rmn4s rmnoinv = 1./rmno ! ! Init fields to be read (not modified by this routine). ! Transform from original f-array to (k,i) arrays: do i=1,len1 nno2k = nno2 nnok = nno nnn2k = nnn2 ntk = nj+nt npnok = nj+npno nps1k = nj+nps nps2k = nj+nps2 nmsk = nj+nms do k=1,kmaxp1 fnno2(k,i) = f(i,nno2k) fnno(k,i) = f(i,nnok ) fnnn2(k,i) = f(i,nnn2k) fnt(k,i) = f(i,ntk ) fnpno(k,i) = f(i,npnok) fnps1(k,i) = f(i,nps1k) fnps2(k,i) = f(i,nps2k) fnms(k,i) = f(i,nmsk) nno2k = nno2k+1 nnok = nnok+1 nnn2k = nnn2k+1 ntk = ntk+1 npnok = npnok+1 nps1k = nps1k+1 nps2k = nps2k+1 nmsk = nmsk+1 beta9ki(k,i) = beta9(i,k) enddo enddo C **** C **** COMPUTE S14 = RSQ, S15 = RSP C **** S1 = TAU(R), S2 = TAU(Q) C **** do i=1,len1 if (idn(i)==1) then do k=1,kmaxp1 s1k(k,i) = SIGEUV(1,49)*fnno2(k,i)+SIGEUV(2,49)*fnno(k,i)+ | SIGEUV(3,49)*fnnn2(k,i) s2k(k,i) = SIGEUV(1,20)*fnno2(k,i)+SIGEUV(2,20)*fnno(k,i)+ | SIGEUV(3,20)*fnnn2(k,i) C **** C **** S3 = TAU(1), S4 = TAU(2), S5 = TAU(3) C **** s3k(k,i) = 1.3*s1k(k,i) s4k(k,i) = 2.0*s1k(k,i) s5k(k,i) = 2.5*s1k(k,i) C **** C **** IF(TAU.GT.9.0) TAU = 9.0 C **** if (s1k(k,i) > 9.) s1k(k,i) = 9. if (s2k(k,i) > 9.) s2k(k,i) = 9. if (s3k(k,i) > 9.) s3k(k,i) = 9. if (s4k(k,i) > 9.) s4k(k,i) = 9. if (s5k(k,i) > 9.) s5k(k,i) = 9. C **** C **** TAU(N) = EXP(-TAU(N)) FOR N = 1,3,1 C **** s3k(k,i) = exp(-s3k(k,i)) s4k(k,i) = exp(-s4k(k,i)) s5k(k,i) = exp(-s5k(k,i)) C **** C **** S6 = EXP(-TAU(R)), S7 = EXP(-TAU(Q)) C **** s6k(k,i) = exp(-s1k(k,i)) s7k(k,i) = exp(-s2k(k,i)) C **** C **** S14 = RSP, S15 = RSQ C **** s14k(k,i) = s6k(k,i)+2.*(s3k(k,i)+s4k(k,i)+s5k(k,i)) s15k(k,i) = 1.5*s6k(k,i)/(s14k(k,i)+s2k(k,i)/s1k(k,i)* | s7k(k,i)) s14k(k,i) = 2.4*s6k(k,i)/s14k(k,i) enddo ! k=1,kmaxp1 endif ! idn==1 enddo ! i=1,len1 C **** C **** S1 = PSI(O2), S2 = PSI(O), S3 = PSI(N2), S4 = PSI(N4S) C **** LEVELS 2 THRU KMAXP1 C **** do i=1,len1 if (idn(i)==1) then do k=1,kmax s1k(k+1,i)=.5*(fnps1(k,i)+fnps1(k+1,i)) s2k(k+1,i)=.5*(fnps2(k,i)+fnps2(k+1,i)) s4k(k+1,i)=0. enddo endif enddo C **** C **** LEVEL 1 bottom boundary C **** do i=1,len1 s1k(1,i) = .5*((b(i,1,1)+1.)*fnps1(1,i)+b(i,1,2)*fnps2(1,i)+ | fb(i,1)) s2k(1,i) = .5*(b(i,2,1)*fnps1(1,i)+(b(i,2,2)+1.)*fnps2(1,i)+ | fb(i,2)) s4k(1,i) = 0. enddo C **** C **** NOW CALCULATE S3 = PSI(N2) C **** do i=1,len1 if (idn(i)==1) then do k=1,kmaxp1 s3k(k,i)=1.-s1k(k,i)-s2k(k,i) enddo endif enddo C **** C **** CALCULATE IONIZATION CONTRIBUTIONS IN S8 - S13 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 **** CONTRIBUTIONS TO Q, J AND QTEF FROM EUV C **** ! Init fields to be summed: do i=1,len1 do k=1,kmaxp1 fnq(k,i) = 0. fnrj(k,i) = 0. fnqtef(k,i) = 0. fnqop2p(k,i) = 0. fnqop2d(k,i) = 0. fnqo2p(k,i) = 0. ! as in old clearnq fnqop(k,i) = 0. ! as in old clearnq fnqn2p(k,i) = 0. ! as in old clearnq fnqnop(k,i) = 0. ! as in old clearnq fnqnp(k,i) = 0. ! as in old clearnq s8k(k,i) = 0. s9k(k,i) = 0. s10k(k,i) = 0. s11k(k,i) = 0. s12k(k,i) = 0. s13k(k,i) = 0. enddo enddo ! ! s[5,6,7]nk must have a lamda dimension because a new sum is needed ! at each wavelength. Note the first dimension is l1:lmax (not lmax). do i=1,len1 if (idn(i)==1) then do k=1,kmaxp1 do n=l1,lmax s5nk(n,k,i) = 0. s6nk(n,k,i) = 0. s7nk(n,k,i) = 0. enddo enddo endif enddo C **** C **** SUMMATION OVER WAVE LENGTH (wavelength loop on inside) C **** 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 **** SUMMATION OVER O2, O, N2 C **** do i=1,zimxp if (idn(i) == 1) then ! daytime do k=1,zkmxp do n=l1,lmax ! 16,59 (44 iters) s5nk(n,k,i)=s5nk(n,k,i) + | sigeuv(1,n)*fnno2(k,i)+ | sigeuv(2,n)*fnno(k,i) + | sigeuv(3,n)*fnnn2(k,i) s6nk(n,k,i)=s6nk(n,k,i)+ | sigeuv(1,n)*s1k(k,i)*rmassinv(1)+ | sigeuv(2,n)*s2k(k,i)*rmassinv(2)+ | sigeuv(3,n)*s3k(k,i)*rmassinv(3) s7nk(n,k,i)=s7nk(n,k,i)+sigmas(1,n)+sigmas(2,n)+ | sigmas(3,n) s5nk(n,k,i) = feuv(n)*exp(-s5nk(n,k,i)) s7nk(n,k,i) = s7nk(n,k,i)+sigmas(4,n) fnq(k,i) = fnq(k,i)+hc*rlmeuvinv(n)*s5nk(n,k,i)* | s6nk(n,k,i) fnrj(k,i) = fnrj(k,i)+s5nk(n,k,i)*s7nk(n,k,i) fnqop2p(k,i) = fnqop2p(k,i)+s5nk(n,k,i)*sigmas(6,n)* | s2k(k,i)*rmassinv(2) fnqop2d(k,i) = fnqop2d(k,i)+s5nk(n,k,i)*sigmas(5,n)* | s2k(k,i)*rmassinv(2) C **** C **** ADD IONIZATION FOR O2, O, N2 TO S8, S9, S10 C **** s8k(k,i) = s8k(k,i)+sigmas(1,n)*s5nk(n,k,i)*s1k(k,i)* | rmassinv(1) s9k(k,i) = s9k(k,i)+sigmas(2,n)*s5nk(n,k,i)*s2k(k,i)* | rmassinv(2) s10k(k,i)=s10k(k,i)+sigmas(3,n)*s5nk(n,k,i)*s3k(k,i)* | rmassinv(3) ! fnqtef(k,i) = fnqtef(k,i)+(sigeuv(3,n)-sigmas(3,n))* | s5nk(n,k,i)*s3k(k,i)*rmassinv(3)*2. ! *2. is a snoe mod s11k(k,i) = s11k(k,i)+sigmas(4,n)*s5nk(n,k,i)*s4k(k,i)* | rmn4sinv s12k(k,i) = s12k(k,i)+sigmas(1,n)*s5nk(n,k,i)*s1k(k,i)* | rmassinv(1)*bro2(n) s13k(k,i) = s13k(k,i)+sigmas(3,n)*s5nk(n,k,i)*s3k(k,i)* | rmassinv(3)*brn2(n) enddo ! n=l1,lmax enddo ! k=1,zkmxp endif ! idn(i)==1 enddo ! i=1,zimxp C **** C **** MULTIPLY Q BY EFFICIENCY FACTOR C **** do i=1,len1 if (idn(i)==1) then do k=1,kmaxp1 fnq(k,i) = fnq(k,i)*euveff(k)*c(85) enddo endif enddo C **** C **** CALCULATE N*MBAR IN S5 C **** S6 = T(TOTAL,K) do i=1,len1 if (idn(i)==1) then s6k(1,i)=fnt(kmaxp1,i)+T0(1) do k=2,kmax s6k(k,i)=.5*(fnt(k-1,i)+fnt(k,i))+T0(K) enddo s6k(kmaxp1,i)=fnt(kmaxp1-1,i)+T0(KMAXP1) endif enddo C **** C **** S7 = EXPS (K) C **** do i=1,len1 if (idn(i) == 1) then ! daytime do k=1,kmax s7k(k,i)=C(87)*EXPS(K) enddo s7k(kmaxp1,i) = c(86)*exps(kmax) endif enddo C **** C **** CALCULATE CONTRIBUTIONS TO NQO2P, NQOP, NQN2P, NQNP C **** AND NQTEF C **** disn2pki = 0. ! whole array init do i=1,len1 if (idn(i) == 1) then ! daytime do k=1,kmaxp1 s5k(k,i)=C(81)*s7k(k,i)*fnms(k,i)/(C(84)*s6k(k,i)) ! N*MBAR disn2pki(k,i) = s10k(k,i)*s14k(k,i)*s5k(k,i) fnqo2p(k,i) = fnqo2p(k,i)+(s8k(k,i)* | (1.+s15k(k,i))-s12k(k,i))*s5k(k,i) fnqn2p(k,i) = fnqn2p(k,i)+(s10k(k,i)* | (1.+s14k(k,i))-s13k(k,i))*s5k(k,i) fnqnp(k,i) = fnqnp(k,i)+(s11k(k,i)+s13k(k,i))*s5k(k,i) fnqtef(k,i) = fnqtef(k,i)*s5k(k,i) s4k(k,i) = fnqop2p(k,i)+fnqop2d(k,i)+s9k(k,i)+s12k(k,i) fnqop2p(k,i) = (fnqop2p(k,i)+s4k(k,i)*s14k(k,i)*0.22)* | s5k(k,i) fnqop2d(k,i) = (fnqop2d(k,i)+s4k(k,i)*s14k(k,i)*0.24)* | s5k(k,i) fnqop(k,i) = fnqop(k,i)+(s9k(k,i)+s12k(k,i)+s4k(k,i)* | s14k(k,i)*0.56)*s5k(k,i) factor=C(85)*C(81)/C(57)*EXPS(1)*C(86)**(2*K-3) s8k(k,i)=factor/((s1k(k,i)*rmassinv(1)+s2k(k,i)*rmassinv(2)+ | s3k(k,i)*rmassinv(3))*s6k(k,i)) s8k(k,i)=s8k(k,i)*(QUENCH(1)*s3k(k,i)*rmassinv(3)+ | QUENCH(2)*s1k(k,i)*rmassinv(1)) s8k(k,i)=QUENCH(3)*s8k(k,i)/(QUENCH(4)+s8k(k,i)) enddo ! k=1,kmaxp1 endif ! idn(i)==1 enddo ! i=1,len1 C **** C **** CALCULATE NO IONIZATION AND ADD TO NQNOP C **** do i=1,len1 if (idn(i) == 1) then ! daytime fnqnop(1,i)=fnqnop(1,i)+beta9ki(1,i)*fnpno(1,i)* | s5k(1,i)*rmnoinv do k=2,kmaxp1 fnqnop(k,i)=fnqnop(k,i)+beta9ki(k,i)*.5*(fnpno(k,i)+ | fnpno(k-1,i))*s5k(k,i)*rmnoinv enddo endif enddo C **** C **** S7=SUM OVER WAVE LENGTH(SIGMA*F*EXP(-SIGMA*CHAPMAN)* C **** (HC/LAMDA-DO2)) C **** INITIALIZE S7 C **** s7k = 0. ! whole array init C **** C **** SUMMATION OVER WAVE LENGTH C **** do i=1,len1 if (idn(i) == 1) then do k=1,kmaxp1 do n=1,l1-1 s9k(k,i)=sigsrc(n)*fsrc(n)*exp(-sigsrc(n)*fnno2(k,i)) s7k(k,i)=s7k(k,i)+s9k(k,i)*(C(60)*rlmsrcinv(n)-do22+ | s8k(k,i)) fnrj(k,i)=fnrj(k,i)+s9k(k,i) enddo fnq(k,i)=fnq(k,i)+s7k(k,i)*C(85)*s1k(k,i)*rmassinv(1) enddo endif enddo C **** C **** UPDATE Q C **** CONTRIBUTIONS FROM SCHUMANN RUNGE BANDS C **** S7=P3*F C **** do i=1,len1 if (idn(i)==1) then do k=1,kmaxp1 s8k(k,i)=SQRT(fnno2(k,i)) if (fnno2(k,i) >= 1.e18) then s7k(k,i) = (1./(aband*fnno2(k,i)+bband*s8k(k,i)))* | (1.+0.11*(f107-65.)/165.)*sfeps else s7k(k,i) = cband*(1.+0.11*(f107-65.)/165.)*sfeps endif fnq(k,i)=fnq(k,i)+s7k(k,i)*C(85)*s1k(k,i)*rmassinv(1)*e3 fnrj(k,i)=fnrj(k,i)+s7k(k,i)/do2 enddo endif enddo ! ! Update original f-array fields: nqk = nq nrjk = nrj nqop2pk = nqop2p nqop2dk = nqop2d nqtefk = nqtef nqo2pk = nqo2p nqopk = nqop nqn2pk = nqn2p nqnpk = nqnp nqnopk = nqnop do k=1,kmaxp1 do i=1,len1 f(i,nqk) = fnq(k,i) f(i,nrjk) = fnrj(k,i) f(i,nqop2pk) = fnqop2p(k,i) f(i,nqop2dk) = fnqop2d(k,i) f(i,nqtefk) = fnqtef(k,i) f(i,nqo2pk) = fnqo2p(k,i) f(i,nqopk) = fnqop(k,i) f(i,nqn2pk) = fnqn2p(k,i) f(i,nqnpk) = fnqnp(k,i) f(i,nqnopk) = fnqnop(k,i) disn2p(i,k) = disn2pki(k,i) enddo ! i=1,len1 nqk = nqk+1 nrjk = nrjk+1 nqop2pk = nqop2pk+1 nqop2dk = nqop2dk+1 nqtefk = nqtefk+1 nqo2pk = nqo2pk+1 nqopk = nqopk+1 nqn2pk = nqn2pk+1 nqnpk = nqnpk+1 nqnopk = nqnopk+1 enddo ! k=1,kmaxp1 return end subroutine qrj !--------------------------------------------------------------- subroutine init_sflux use input_module,only: f107,f107a use init_module,only: sfeps ! ! Flux initialization once per time step, called from advnce. ! (this call moved from qrj in multi-task version) ! include "params.h" include "cons.h" include "qrj_coeff.h" real :: EUVFLX(37),wave1(lmax),wave2(lmax) ! output from euvac real :: flya,hlybr,fexvir,hlya,heiew,xuvfac integer :: iscale,n ! ! 2/00: ! Use f107,f107a from input_mod (these are either provided by ! the user, or obtained from GPI database, see input_mod.f and ! gpi_mod.f). ! HLYBR = 0. FEXVIR = 0. HLYA = 3.E+11+0.4E+10*(f107-70.) HEIEW = 0. C XUVFAC =0. C XUVFAC = 2.0 - (f107-68.0) / (243.0-68.0) XUVFAC = 4.0 - (f107-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)*SFEPS 4 CONTINUE DO 5 N = 1,15 FSRC(N) = SFLUX(N)*SFEPS 5 CONTINUE DO 55 N=1,37 NN = N+15 FEUV(NN) = EUVFLX(N)*SFEPS 55 CONTINUE return end !--------------------------------------------------------------- block data set_qrj_coeff include "params.h" include "qrj_coeff.h" 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 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" include "qrj_coeff.h" 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) ! 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