SUBROUTINE HYDLIW(ID,ABSOH,EMISH) C ================================= C C opacity and emissivity of hydrogen lines C INCLUDE 'PARAMS.FOR' INCLUDE 'MODELP.FOR' INCLUDE 'SYNTHP.FOR' INCLUDE 'WINCOM.FOR' PARAMETER (FRH1=3.28805E15,FRH2=FRH1/4.,UN=1.,SIXTH=1./6.) PARAMETER (CPP=4.1412E-16,CPJ=157803.) PARAMETER (C00=1.25E-9,CDOP=1.284523E12,CID=0.02654,TWO=2.) PARAMETER (CPJ4=CPJ/4.,AL10=2.3025851,CINV=UN/2.997925E18) PARAMETER (CID1=0.01497) common/lasers/lasdel common/quasun/nunalp,nunbet,nungam,nunbal DIMENSION PJ(40),PRF0(54),OSCH(4,22), * ABSO(MFREQ),EMIS(MFREQ),ABSOH(MFREQ),EMISH(MFREQ) DATA FRH /3.289017E15/ DATA INIT /0/ C if(iath.le.0) return izz=1 C IF(INIT.EQ.0) THEN DO I=1,4 DO J=I+1,22 CALL STARK0(I,J,IZZ,XK,WL0,FIJ,FIJ0) WLINE(I,J)=WL0 OSCH(I,J)=FIJ+FIJ0 END DO END DO INIT=1 END IF DO IJ=1,NFREQ ABSO(IJ)=0. EMIS(IJ)=0. ABSOH(IJ)=0. EMISH(IJ)=0. END DO T=TEMP(ID) T1=UN/T SQT=SQRT(T) ANE=ELEC(ID) ANES=EXP(SIXTH*LOG(ANE)) C C populations of the first 40 levels of hydrogen C ANP=POPUL(NKH,ID) PP=CPP*ANE*ANP*T1/SQT NLH=N1H-N0HN+1 if(ifwop(n1h).lt.0) nlh=nlh-1 DO 5 IL=1,40 X=IL*IL IF(IL.LE.NLH) PJ(IL)=POPUL(N0HN+IL-1,ID) IF(IL.GT.NLH) PJ(IL)=PP*EXP(CPJ/X*T1)*X*wnhint(il,id) 5 CONTINUE p2=pp*exp(cpj4*t1)*4.*wnhint(2,id) C C Frequency- and line-independent parameters for evaluating the C asymptotic Stark profile C F00=C00*ANES*ANES*ANES*ANES DOP0=1.E8*SQRT(1.65E8*T+VTURB(ID)) C C ------------------------------------------------------------------- C overall loop over spectral series (only in the infrared region) C ------------------------------------------------------------------- C DO 300 IJ=1,NFREQ IF(IHYLW(IJ).LE.0) GO TO 300 ISERL=ILOWHW(IJ) ISERU=ILOWHW(IJ) IF(WLAM(IJ).GT.17000..AND.WLAM(IJ).LE.21000.) THEN ISERL=3 ISERU=4 ELSE IF(WLAM(IJ).GT.22700..AND.WLAM(IJ).LE.29000.) THEN ISERL=4 ISERU=5 ELSE IF(WLAM(IJ).GT.32800..AND.WLAM(IJ).LE.37000.) THEN ISERL=5 ISERU=6 ELSE IF(WLAM(IJ).GT.37000..AND.WLAM(IJ).LE.44600.) THEN ISERL=4 ISERU=6 ELSE IF(WLAM(IJ).GT.44660..AND.WLAM(IJ).LE.58300.) THEN ISERL=5 ISERU=7 ELSE IF(WLAM(IJ).GT.58300..AND.WLAM(IJ).LE.72000.) THEN ISERL=6 ISERU=8 ELSE IF(WLAM(IJ).GT.72000..AND.WLAM(IJ).LE.73800.) THEN ISERL=5 ISERU=8 ELSE IF(WLAM(IJ).GT.73800..AND.WLAM(IJ).LE.77000.) THEN ISERL=5 ISERU=9 ELSE IF(WLAM(IJ).GT.77000.) THEN ISERL=6 ISERU=9 END IF C if(iserl.eq.3.and.iseru.eq.3.and.nunbal.gt.0) iserl=2 C ABSO(IJ)=0. EMIS(IJ)=0. DO 200 I=ISERL,ISERU II=I*I XII=UN/II PLTEI=PP*EXP(CPJ*T1*XII)*II POPI=PJ(I) IF(I.EQ.1) FRH=3.28805E15 C C determination of which hydrogen lines contribute in a current C frequency region C M1=M10W(IJ) IF(I.LT.ILOWHW(IJ)) M1=ILOWHW(IJ)-1 M2=M1+1 IF(M1.LT.I+1) M1=I+1 IF(grav.lt.3..and.M1.LE.16.AND.I.EQ.7) GO TO 10 IF(grav.lt.3..and.M1.LE.14.AND.I.EQ.6) GO TO 10 IF(grav.lt.3..and.M1.LE.12.AND.I.EQ.5) GO TO 10 IF(grav.lt.3..and.M1.LE.10.AND.I.EQ.4) GO TO 10 IF(grav.lt.3..and.M1.LE.8.AND.I.EQ.3) GO TO 10 IF(grav.lt.3..and.M1.LE.6.AND.I.EQ.2) GO TO 10 IF(grav.lt.3..and.M1.LE.4.AND.I.EQ.1) GO TO 10 M1=M1-1 M2=M20W(IJ)+3 IF(M1.LT.I+1) M1=I+1 10 CONTINUE if(grav.gt.3.) then m2=m2+5 m1=m1-3 if(m1.gt.i+6) m1=m1-3 end if if(grav.gt.6.) then m2=m2+2 m1=m1-1 if(m1.gt.i+6) m1=m1-1 end if IF(M1.LT.I+1) M1=I+1 c if(m2.gt.30) then c m2=m20W(IJ)+8 c m1=m1-4 c end if IF(M2.GT.40) M2=40 c if(id.eq.1) write(6,666) i,m1,m2 c 666 format(/' hydrogen lines contribute - ilow=',i2,', iup from ',i3, c * ' to',i3/) C A=0. E=0. C C loop over lines which contribute at given wavelength region C DO 100 J=M1,M2 IF(I.EQ.1.AND.J.LE.5.AND.IOPHLI.LT.0) GO TO 100 ILINE=0 JJ=J*J XJJ=UN/JJ ABTRA=PJ(I)*WNHINT(J,ID) EMTRA=PJ(J)*WNHINT(I,ID)*II*XJJ*EXP(CPJ*(XII-XJJ)*T1) if(i.le.2.and.j.le.i+2) then abtra=pj(i) emtra=pj(j)*wnhint(i,id)/wnhint(j,id)* * ii*xjj*exp(cpj*(xii-xjj)*t1) end if IF(I.LE.4.AND.J.LE.22) ILINE=ILIN0(I,J) c c quasi-molecular opacity for Lyman-alpha and beta satellites c lquasi=i.eq.1.and.j.eq.2.and.nunalp.gt.0 lquasi=lquasi.or.i.eq.1.and.j.eq.3.and.nunbet.gt.0 lquasi=lquasi.or.i.eq.1.and.j.eq.4.and.nungam.gt.0 lquasi=lquasi.or.i.eq.2.and.j.eq.3.and.nunbal.gt.0 if(lquasi) then CALL STARK0(I,J,izz,XKIJ,WL0,FIJ,FIJ0) FXK=F00*XKIJ FXK1=UN/FXK DOP=DOP0/WL0 DBETA=WL0*WL0*CINV*FXK1 BETAD=DOP*DBETA FID=CID*FIJ*DBETA CALL DIVSTR(AD,DIV) fr=freq(ij) BETA=ABS(WLAM(IJ)-WL0)*FXK1 call allard(wlam(ij),popi,anp,sg,i,j) sg=sg+STARKA(BETA,AD,DIV,UN)*FID ABSO(IJ)=ABSO(IJ)+SG*ABTRA EMIS(IJ)=EMIS(IJ)+SG*EMTRA go to 100 end if c c lines with special Stark broadening tables c IF(ILINE.GT.0) THEN NWL=NWLHYD(ILINE) DO IWL=1,NWL PRF0(IWL)=PRFHYD(ILINE,ID,IWL) END DO FID=CID*OSCH(I,J) AL=ABS(WLAM(IJ)-WLINE(I,J)) IF(AL.LT.1.E-4) AL=1.E-4 IF(ILEMKE.EQ.1) AL=AL/F00 AL=LOG10(AL) DO 30 IWL=1,NWL-1 IW0=IWL IF(AL.LE.WLHYD(ILINE,IWL+1)) GO TO 40 30 CONTINUE 40 IW1=IW0+1 PRFF=(PRF0(IW0)*(WLHYD(ILINE,IW1)-AL)+PRF0(IW1)* * (AL-WLHYD(ILINE,IW0)))/ * (WLHYD(ILINE,IW1)-WLHYD(ILINE,IW0)) SG=EXP(PRFF*AL10)*FID IF(ILEMKE.EQ.1) SG=SG*WLINE(I,J)**2*CINV/F00 ABSO(IJ)=ABSO(IJ)+SG*ABTRA EMIS(IJ)=EMIS(IJ)+SG*EMTRA c c lines without special Stark broadening tables c ELSE CALL STARK0(I,J,izz,XKIJ,WL0,FIJ,FIJ0) FXK=F00*XKIJ FXK1=UN/FXK DOP=DOP0/WL0 DBETA=WL0*WL0*CINV*FXK1 BETAD=DOP*DBETA FID=CID*FIJ*DBETA CALL DIVSTR(AD,DIV) fr=freq(ij) BETA=ABS(WLAM(IJ)-WL0)*FXK1 SG=STARKA(BETA,AD,DIV,TWO)*FID if(iophli.eq.2.and.i.eq.1.and.j.eq.2) * sg=sg*feautr(fr,id) ABSO(IJ)=ABSO(IJ)+SG*ABTRA EMIS(IJ)=EMIS(IJ)+SG*EMTRA END IF 100 CONTINUE 200 CONTINUE C C ---------------------------- C total opacity and emissivity C ---------------------------- C F=FREQ(IJ) F15=F*1.E-15 XKF=EXP(-4.79928e-11*F*T1) XKFB=XKF*1.4743E-2*F15*F15*F15 if(abso(ij).le.0. .and. lasdel) then abso(ij)=0. emis(ij)=0. endif ABSOH(IJ)=ABSO(IJ)-XKF*EMIS(IJ) EMISH(IJ)=XKFB*EMIS(IJ) 300 CONTINUE RETURN END