diff --git a/cluster/edfb.f b/cluster/edfb.f
new file mode 100644
index 0000000000000000000000000000000000000000..e54ab7462955ef05e63f11669a2357ca7edffec2
--- /dev/null
+++ b/cluster/edfb.f
@@ -0,0 +1,336 @@
+ PROGRAM EDFB
+CCC 160630
+CCC IES=1 FOR SURROUNDING EXTERNAL SPHERE CENTERED AT ORIGIN;
+CCC HOMOGENEOUS MATERIAL WITHIN THE EXT. SPHERE IS READ IN
+CCC AS THAT OF THE (NSHL+1)-TH LAYER OF THE 1-ST SPHERE:
+CCC RCF(1,NSHL)=1.0D0, RCF(1,NSHL+1)>1.0D0 AND
+CCC ROS(1) IS RADIUS OF SPHERE 1
+CCC
+CCC READ DATA FOR BUILDING VECTOR XIV FROM WITHIN SUB INXI
+CCC
+CCC IDFC>0 WHEN ALL DIEL. FUNCT. OF SPHERES ARE CONSTANTS;
+CCC IDFC=0 WHEN DIEL. FUNCT. OF SPHERES DEPEND ON XI;
+CCC IDFC<0 WHEN DIEL. FUNCT. OF SPHERES ARE AT XIP VALUE ONLY AND
+CCC XI IS SCALE FACTOR FOR DIMENSIONS
+CCC INSN CHOOSES THE VARIABLE THE DIEL. FUNCT. DEPEND ON:
+CCC (INSN=1)=XI;
+CCC (INSN=2)=WN (THE WAVENUMBER, IN m**-1);
+CCC (INSN=3)=WL (THE WAVELENGTH, IN m);
+CCC (INSN=4)=PU (THE ANGULAR FREQUENCY, IN s**-1);
+CCC (INSN=5)=EV (THE PHOTON ENERGY, IN ev);
+CCC INSTPC>0 WHEN VARIABLE INCREASES WITH A CONSTANT STEP;
+CCC INSTPC=0 WHEN VARIABLE IS SAVED IN A VECTOR YOU READ BEFORE
+CCC SUPPORTS EXPERIMENTAL DIELECTRIC FUNCTIONS ONLY
+CCC NSPH=6
+ IMPLICIT REAL*8(A-H,O-Z)
+CCC COMMON/C1/DC0(NSHL-NTL+1),DC0M(NSHL-NTL+1,NSPH,NXI),
+CCC 1ROS(NSPH),RCF(NSPH,NSHL+1),IOG(NSPH),NSHL(NSPH)
+ COMMON/C1/DC0(5),DC0M(5,6,200),
+ 1ROS(6),RCF(6,9),IOG(6),NSHL(6)
+ COMPLEX*16 DC0,DC0M
+CCC COMMON/C3/XIV(NXI),WNS(NXI),WLS(NXI),PUS(NXI),EVS(NXI),
+CCC 1VSS(NXI),VNS(5)
+ COMMON/C3/XIV(200),WNS(200),WLS(200),PUS(200),EVS(200),
+ 1VSS(200),VNS(5)
+ CHARACTER*3 VNS
+ 5010 FORMAT(16I5)
+ 6005 FORMAT(' SPHERE N.',I4)
+ 6009 FORMAT(' NONTRANSITION LAYER N.',I2,', SCALE = ',A3)
+ 6010 FORMAT(I5,1X,1PD12.4,1PD12.4)
+ IR=5
+ IW=6
+ IT=7
+CCC
+CCC SIZE(I)=VK*ROS(I) IS SIZE PARAMETER in vacuo
+CCC
+CCC
+CCC READING OF DIELECTRIC FUNCTIONS DRIVEN BY ICI DEFINED BELOW
+CCC
+ OPEN(IR,FILE='DEDFB',STATUS='OLD')
+ READ(IR,*)NSPH,IES
+ IF(IES.NE.0)IES=1
+ READ(IR,*)EXDC,WP,XIP,IDFC,NXI,INSTPC,INSN
+ OPEN(IW,FILE='OEDFB',STATUS='UNKNOWN')
+ CALL INXI(IR,IW,WP,XIP,IDFC,NXI,INSTPC,INSN)
+ READ(IR,5010)(IOG(I),I=1,NSPH)
+ DO 113 I=1,NSPH
+ IF(IOG(I).LT.I)GO TO 113
+ READ(IR,*)NSHL(I),ROS(I)
+ NSH=NSHL(I)
+ IF(I.EQ.1)NSH=NSH+IES
+ DO 112 NS=1,NSH
+ 112 READ(IR,*)RCF(I,NS)
+ 113 CONTINUE
+ OPEN(IT,FILE='TEDF',FORM='UNFORMATTED',STATUS='UNKNOWN')
+ WRITE(IT)NSPH
+ WRITE(IT)(IOG(I),I=1,NSPH)
+ WRITE(IT)EXDC,WP,XIP,IDFC,NXI
+ WRITE(IT)(XIV(I),I=1,NXI)
+ DO 115 I=1,NSPH
+ IF(IOG(I).LT.I)GO TO 115
+ WRITE(IT)NSHL(I),ROS(I)
+ NSH=NSHL(I)
+ IF(I.EQ.1)NSH=NSH+IES
+ WRITE(IT)(RCF(I,NS),NS=1,NSH)
+ 115 CONTINUE
+ DO 468 JXI=1,NXI
+ IF((IDFC.NE.0).AND.(JXI.GT.1))GO TO 468
+ DO 162 I=1,NSPH
+ IF(IOG(I).LT.I)GO TO 162
+CCC
+ NSH=NSHL(I)
+ ICI=(NSH+1)/2
+ IF(I.EQ.1)ICI=ICI+IES
+ DO 157 IC=1,ICI
+ READ(IR,*)DC0(IC)
+ 157 DC0M(IC,I,JXI)=DC0(IC)
+CCC
+ WRITE(IT)(DC0(IC),IC=1,ICI)
+ 162 CONTINUE
+ 468 CONTINUE
+ IF(IDFC.EQ.0)GO TO 474
+ WRITE(IW,*)' DIELECTRIC CONSTANTS'
+ DO 473 I=1,NSPH
+ IF(IOG(I).NE.I)GO TO 473
+ ICI=(NSHL(I)+1)/2
+ IF(I.EQ.1)ICI=ICI+IES
+ WRITE(IW,6005)I
+ DO 472 IC=1,ICI
+ WRITE(IW,6010)IC,DC0M(IC,I,1)
+ 472 CONTINUE
+ 473 CONTINUE
+ GO TO 499
+ 474 WRITE(IW,*)' DIELECTRIC FUNCTIONS'
+ DO 478 I=1,NSPH
+ IF(IOG(I).NE.I)GO TO 478
+ ICI=(NSHL(I)+1)/2
+ IF(I.EQ.1)ICI=ICI+IES
+ WRITE(IW,6005)I
+ DO 477 IC=1,ICI
+ WRITE(IW,6009)IC,VNS(INSN)
+ DO 476 JXI=1,NXI
+ WRITE(IW,6010)JXI,DC0M(IC,I,JXI)
+ 476 CONTINUE
+ 477 CONTINUE
+ 478 CONTINUE
+ 499 CLOSE(IR)
+ CLOSE(IW)
+ CLOSE(IT)
+ STOP
+ END
+ SUBROUTINE INXI(IR,IW,WP,XIP,IDFC,NXI,INSTPC,INSN)
+ IMPLICIT REAL*8(A-H,O-Z)
+ COMMON/C3/XIV(200),WNS(200),WLS(200),PUS(200),EVS(200),
+ 1VSS(200),VNS(5)
+ CHARACTER*3 VNS
+ PIGT=DACOS(0.0D0)*4.0D0
+ EVC=6.5821188D-16
+ IF(IDFC.LT.0)GO TO 300
+ IF(INSTPC.EQ.0)GO TO 200
+CCC VLST=V+(NXI-1)*VSTP
+ GO TO(105,125,145,165,185),INSN
+ RETURN
+ 105 READ(IR,*)XI,XISTP
+ DO 110 JXI=1,NXI
+ PU=XI*WP
+ WN=PU/3.0D08
+ VNS(INSN)='XIV'
+ VSS(JXI)=XI
+ XIV(JXI)=XI
+ PUS(JXI)=PU
+ EVS(JXI)=PU*EVC
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ 110 XI=XI+XISTP
+ WRITE(IW,6601)
+ 6601 FORMAT
+ 1(2X,'JXI',5X,'XIV',10X,'WNS',10X,'WLS',10X,'PUS',10X,'EVS')
+ WRITE(IW,6600)
+ 1(JXI,XIV(JXI),WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),JXI=1,NXI)
+ 6600 FORMAT((I5,5(1PD13.4)))
+ RETURN
+ 125 READ(IR,*)WN,WNSTP
+ DO 130 JXI=1,NXI
+ XI=3.0D08*WN/WP
+ PU=XI*WP
+ VNS(INSN)='WNS'
+ VSS(JXI)=WN
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ XIV(JXI)=XI
+ PUS(JXI)=PU
+ EVS(JXI)=PU*EVC
+ 130 WN=WN+WNSTP
+ WRITE(IW,6602)
+ 6602 FORMAT
+ 1(2X,'JXI',5X,'WNS',10X,'WLS',10X,'PUS',10X,'EVS',10X,'XIV')
+ WRITE(IW,6600)
+ 1(JXI,WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 145 READ(IR,*)WL,WLSTP
+ DO 150 JXI=1,NXI
+ WN=PIGT/WL
+ XI=3.0D08*WN/WP
+ PU=XI*WP
+ VNS(INSN)='WLS'
+ VSS(JXI)=WL
+ WLS(JXI)=WL
+ WNS(JXI)=WN
+ XIV(JXI)=XI
+ PUS(JXI)=PU
+ EVS(JXI)=PU*EVC
+ 150 WL=WL+WLSTP
+ WRITE(IW,6603)
+ 6603 FORMAT
+ 1(2X,'JXI',5X,'WLS',10X,'WNS',10X,'PUS',10X,'EVS',10X,'XIV')
+ WRITE(IW,6600)
+ 1(JXI,WLS(JXI),WNS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 165 READ(IR,*)PU,PUSTP
+ DO 170 JXI=1,NXI
+ XI=PU/WP
+ WN=PU/3.0D08
+ VNS(INSN)='PUS'
+ VSS(JXI)=PU
+ PUS(JXI)=PU
+ XIV(JXI)=XI
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ EVS(JXI)=PU*EVC
+ 170 PU=PU+PUSTP
+ WRITE(IW,6604)
+ 6604 FORMAT
+ 1(2X,'JXI',5X,'PUS',10X,'WNS',10X,'WLS',10X,'EVS',10X,'XIV')
+ WRITE(IW,6600)
+ 1(JXI,PUS(JXI),WNS(JXI),WLS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 185 READ(IR,*)EV,EVSTP
+ DO 190 JXI=1,NXI
+ PU=EV/EVC
+ XI=PU/WP
+ WN=PU/3.0D08
+ VNS(INSN)='EVS'
+ VSS(JXI)=EV
+ EVS(JXI)=EV
+ PUS(JXI)=PU
+ XIV(JXI)=XI
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ 190 EV=EV+EVSTP
+ WRITE(IW,6605)
+ 6605 FORMAT
+ 1(2X,'JXI',5X,'EVS',10X,'WNS',10X,'WLS',10X,'PUS',10X,'XIV')
+ WRITE(IW,6600)
+ 1(JXI,EVS(JXI),WNS(JXI),WLS(JXI),PUS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 200 GO TO(205,225,245,265,285),INSN
+ RETURN
+ 205 DO 210 JXI=1,NXI
+ READ(IR,*)XI
+ PU=XI*WP
+ WN=PU/3.0D08
+ VNS(INSN)='XIV'
+ VSS(JXI)=XI
+ XIV(JXI)=XI
+ PUS(JXI)=PU
+ EVS(JXI)=PU*EVC
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ 210 CONTINUE
+ WRITE(IW,6601)
+ WRITE(IW,6600)
+ 1(JXI,XIV(JXI),WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),JXI=1,NXI)
+ RETURN
+ 225 DO 230 JXI=1,NXI
+ READ(IR,*)WN
+ XI=3.0D08*WN/WP
+ PU=XI*WP
+ VNS(INSN)='WNS'
+ VSS(JXI)=WN
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ XIV(JXI)=XI
+ PUS(JXI)=PU
+ EVS(JXI)=PU*EVC
+ 230 CONTINUE
+ WRITE(IW,6602)
+ WRITE(IW,6600)
+ 1(JXI,WNS(JXI),WLS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 245 DO 250 JXI=1,NXI
+ READ(IR,*)WL
+ WN=PIGT/WL
+ XI=3.0D08*WN/WP
+ PU=XI*WP
+ VNS(INSN)='WLS'
+ VSS(JXI)=WL
+ WLS(JXI)=WL
+ WNS(JXI)=WN
+ XIV(JXI)=XI
+ PUS(JXI)=PU
+ EVS(JXI)=PU*EVC
+ 250 CONTINUE
+ WRITE(IW,6603)
+ WRITE(IW,6600)
+ 1(JXI,WLS(JXI),WNS(JXI),PUS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 265 DO 270 JXI=1,NXI
+ READ(IR,*)PU
+ XI=PU/WP
+ WN=PU/3.0D08
+ VNS(INSN)='PUS'
+ VSS(JXI)=PU
+ PUS(JXI)=PU
+ XIV(JXI)=XI
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ EVS(JXI)=PU*EVC
+ 270 CONTINUE
+ WRITE(IW,6604)
+ WRITE(IW,6600)
+ 1(JXI,PUS(JXI),WNS(JXI),WLS(JXI),EVS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 285 DO 290 JXI=1,NXI
+ READ(IR,*)EV
+ PU=EV/EVC
+ XI=PU/WP
+ WN=PU/3.0D08
+ VNS(INSN)='EVS'
+ VSS(JXI)=EV
+ EVS(JXI)=EV
+ PUS(JXI)=PU
+ XIV(JXI)=XI
+ WNS(JXI)=WN
+ WLS(JXI)=PIGT/WN
+ 290 CONTINUE
+ WRITE(IW,6605)
+ WRITE(IW,6600)
+ 1(JXI,EVS(JXI),WNS(JXI),WLS(JXI),PUS(JXI),XIV(JXI),JXI=1,NXI)
+ RETURN
+ 300 IF(INSTPC.GT.0)GO TO 315
+ DO 310 JXI=1,NXI
+ READ(IR,*)XI
+ VNS(INSN)='XIV'
+ VSS(JXI)=XI
+ XIV(JXI)=XI
+ 310 CONTINUE
+ GO TO 330
+ 315 READ(IR,*)XI,XISTP
+ DO 320 JXI=1,NXI
+ VNS(INSN)='XIV'
+ VSS(JXI)=XI
+ XIV(JXI)=XI
+ 320 XI=XI+XISTP
+ 330 PU=XIP*WP
+ WN=PU/3.0D08
+ WRITE(IW,6611)
+ 6611 FORMAT
+ 1(10X,'XIP',10X,'WN ',10X,'WL ',10X,'PU ',10X,'EV')
+ WRITE(IW,6610)XIP,WN,PIGT/WN,PU,PU*EVC
+ 6610 FORMAT((5X,5(1PD13.4)))
+ WRITE(IW,*)' SCALE FACTORS XI'
+ WRITE(IW,6612)(JXI,XIV(JXI),JXI=1,NXI)
+ 6612 FORMAT(I5,1PD13.4)
+ RETURN
+ END
+CCC
\ No newline at end of file