this is my program below:

SUBROUTINE DFLUX(FLUX,SOL,KSTEP,KINC,TIME,NOEL,NPT,COORDS,

1 JLTYP,TEMP,PRESS,SNAME)

C

INCLUDE 'ABA_PARAM.INC'

C

DIMENSION FLUX(2), TIME(2), COORDS(3)

CHARACTER*80 SNAME

CHARACTER*80 NAME01, NAME02, NAME03

NAME01 = 'ASSEMBLY_Part-2-1_Surf-notrh'

NAME02 = 'ASSEMBLY_Part-2-1_Surf-south'

NAME03 = 'ASSEMBLY_Part-2-1_Surf-top'

NAME04 = 'ASSEMBLY_Part-2-1_Surf-ground'

IF TIME(1) .LE. 24 THEN

REAL TAO = ((12-TIME(1))*15)*3.14/180

ELSE IF (TIME(1) .LE. 48).AND.(TIME(1) .GT. 24) THEN

REAL TAO = ((12-(TIME(1)-24))*15)*3.14/180

ELSE IF (TIME(1) .LE. 72).AND.(TIME(1) .GT. 48) THEN

REAL TAO = ((12-(TIME(1)-48))*15)*3.14/180

ELSE IF (TIME(1) .LE. 96).AND.(TIME(1) .GT. 72) THEN

REAL TAO = ((12-(TIME(1)-72))*15)*3.14/180

ELSE IF (TIME(1) .LE. 120).AND.(TIME(1) .GT. 96) THEN

REAL TAO = ((12-(TIME(1)-96))*15)*3.14/180

C 太阳时角弧度(选取计算日期为2020年7月23日~7月27日，共计5天时间）

REAL Idn=(1325.231*0.66**(1/(0.723*COS(TAO)+0.217))) *

1(0.281-0.525*cos(TAO))

REAL Ids=(1325.231*0.66**(1/(0.723*COS(TAO)+0.217))) *

1(-0.237+0.674*cos(TAO))

REAL Idh=1325.231*0.66**(1/(0.723*COS(TAO)+0.217))

C 以上为太阳直接辐射(w/m2)

REAL Idhh=(359.14-389.62*0.66**(1/(0.723*COS(TAO)+0.217)))*

1(0.723*COS(TAO)+0.217)

REAL Idhn=0.551*Idhh

REAL Idhs=0.551*Idhh

C 以上为考虑太阳散射辐射(w/m2)

REAL Irbn=0.0898*(1325.231*0.66**(1/(0.723*COS(TAO)+0.217))*

1(0.723*COS(TAO)+0.217)+Idhh)

REAL Irbs=Irbn

C 以上为考虑太阳地面反射(w/m2)

IF (SNAME.EQ.NAME01) THEN

FLUX(1)=0.5*(Idn+Idhn+Irbn)*3600

ELSE IF (SNAME.EQ.NAME02) THEN

FLUX(1)=0.5*(Ids+Idhs+Irbs)*3600

ELSE IF (SNAME.EQ.NAME03) THEN

FLUX(1)=0.5*(Idh+Idhh)*3600

ELSE IF (SNAME.EQ.NAME03) THEN

FLUX(1)=0.5*(Idh+Idhh)*3600

ELSE

FLUX(1)=1325

END IF

RETURN

END