Upload New File

sphere/edfb.f

+      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