SpectraRust/tlusty/extracted/corrwm.f
2026-03-19 14:05:33 +08:00

126 lines
3.5 KiB
Fortran

SUBROUTINE CORRWM
C =================
C
C The routine for management of various flags for treating
C frequency points; in particular those connected to the so-called
C "subtraction weights" (in the non-overlapping mode only)
C
INCLUDE 'IMPLIC.FOR'
INCLUDE 'BASICS.FOR'
INCLUDE 'ATOMIC.FOR'
INCLUDE 'MODELQ.FOR'
PARAMETER (T15=1.D-15)
C
NFREQE=0
DO 10 IJ=1,NFREQ
IJEX(IJ)=0
DO ID=1,ND
LSKIP(ID,IJ)=.FALSE.
END DO
if(ifprad.eq.0) then
do id=1,nd
lskip(id,ij)=.true.
end do
end if
IF(IJALI(IJ).NE.0) GO TO 10
NFREQE=NFREQE+1
IJEX(IJ)=NFREQE
IJFR(NFREQE)=IJ
10 CONTINUE
c
if(ifryb.ne.0) then
nfreqe=0
do ij=1,nfreq
ijex(ij)=0
end do
end if
C
IF(IBFINT.LE.0) THEN
DO IJ=1,NFREQ
IJBF(IJ)=IJ
AIJBF(IJ)=UN
END DO
ELSE
IF(ISPODF.EQ.0) THEN
DO IJ=1,NFREQC
IJBF(IJ)=IJ
AIJBF(IJ)=UN
END DO
IF(NFREQ.GT.NFREQC) THEN
DO IJ=NFREQC+1,NFREQ
FR=FREQ(IJ)
IJ0=1
DO IJT=1,NFREQC
IF(FREQ(IJT).LE.FR) THEN
IJ0=IJT
GO TO 12
END IF
END DO
12 IJ1=IJ0-1
A1=(FR-FREQ(IJ0))/(FREQ(IJ1)-FREQ(IJ0))
IJBF(IJ)=IJ1
AIJBF(IJ)=A1
END DO
END IF
ELSE
DO IJ=1,NFREQC-1
IJ0=IFREQB(IJ)
IJ1=IFREQB(IJ+1)
DO KJ=IJ0,IJ1-1
IJBF(KJ)=IJ
AIJBF(KJ)=(FREQ(KJ)-FREQ(IJ1))/(FREQ(IJ0)-FREQ(IJ1))
END DO
END DO
IJ0=IFREQB(NFREQC)
IJBF(IJ0)=NFREQC
AIJBF(IJ0)=UN
END IF
END IF
C
if(nfreqe.gt.mfrex) CALL QUIT('nfreqe.gt.mfrex',nfreqe,mfrex)
C
DO 100 ITR=1,NTRANS
IF(.NOT.LINE(ITR)) GO TO 100
C
C first set up array LSKIP(ID,IJ), which has values
C TRUE - if the radiation at frequency point IJ does not contribute
C radiation pressure (ie. this point belongs to a transition
C for which the user required the radiation pressure to be
C skipped - IABS(INDEXP) chosen as 9 or 19)
C FALSE - normal calculatitn of radiation pressure
C
INX=IABS(INDEXP(ITR))
IF(INX.EQ.9.OR.INX.GE.19) THEN
DO IJ=IFR0(ITR),IFR1(ITR)
DO ID=1,ND
LSKIP(ID,IJ)=.TRUE.
END DO
END DO
END IF
100 CONTINUE
C
IF(NFREQE.GT.0) WRITE(6,609)
DO 110 IJ=1,NFREQ
FR15=FREQ(IJ)*T15
W0E(IJ)=W(IJ)*PI4H/FREQ(IJ)
BNUE(IJ)=BN*FR15*FR15*FR15
IF(IJALI(IJ).NE.0.or.ifryb.gt.0) GO TO 110
if(ispodf.eq.0) then
WRITE(6,610) IJ,FREQ(IJ),W(IJ),PROF(IJ)
else
WRITE(6,610) IJ,FREQ(IJ),W(IJ)
end if
110 CONTINUE
C
DO IJ=1,NFREQ
WC(IJ)=W(IJ)
IF(IJALI(IJ).LE.0) WC(IJ)=0.
END DO
C
609 FORMAT(1H0//' FREQUENCY POINTS AND WEIGHTS - EXPLICIT'/
* ' ---------------------------------------'//
* ' IJ',7X,'FREQ',13X,'WEIGHT',11X,'PROF'/)
610 FORMAT(1H ,I8,1P2D17.8,D15.5,D17.8)
RETURN
END