📄
字号:
program TM2
! -----------------------------------------------------------------
! -----------------------------------------------------------------
! 矩形格子2维光子晶体的时域有限差分FORTRAN程序
!
! -----------------------------------------------------------------
! -----------------------------------------------------------------
!implicit none
external even_jiedian
!***********************************************************
integer i,j,j0,jw
integer nx,ny
integer nt,n
integer nx_colomn,ny_colomn
integer ios
integer grid_vs_colomn2
integer temp1,temp2
integer temp11,temp22
!integer tempp
real*8 even_jiedian
real*8 dx,dy
real*8 dt
real*8 zi_jie_shu
real*8 w,fff
real*8 t
!real*8 maxx
integer nsx1,nsx2
integer nsy1,nsy2
integer n0,m0,iii
character*5 wave
!character*27 filename
!character*31 filename1
character*10 numbers
real*8 jiedian1
real*8 jiedian2
real(8) length_crystal
real(8) width_crystal
real*8 v_b
real*8 v_a
real*8 r
real*8 r_vs_a
real*8 r_vs_b ,light_speed
! 角点
real*8 dd,sina1,cosa1,frad
!角点完成
!能量是否守恒
real*8 sru,schu,szhong,schu0,sru0
real*8 schu01,sru01
!能量守恒定义部分完成
real*8, dimension(:,:),allocatable::u
real*8, dimension(:,:),allocatable::v
real*8, dimension(:,:),allocatable::ez1
real*8, dimension(:,:),allocatable::ez2
real*8, dimension(:,:),allocatable::ez
real*8, dimension(:,:),allocatable::hx
real*8, dimension(:,:),allocatable::hy
real*8, dimension(:,:),allocatable::jiedian
real*8, dimension(:,:),allocatable::jiedian11
real*8, dimension(: ),allocatable::xt
real*8, dimension(: ),allocatable::yt
real*8, dimension(: ),allocatable::xt1
real*8, dimension(: ),allocatable::yt1
real*8, dimension(: ),allocatable::xt2
real*8, dimension(: ),allocatable::yt2
real*8, dimension(: ),allocatable::xt3
real*8, dimension(: ),allocatable::yt3
real*8, dimension(: ),allocatable::xt4
real*8, dimension(: ),allocatable::yt4
real*8, dimension(: ),allocatable::xt5
real*8, dimension(: ),allocatable::yt5
real*8, dimension(: ),allocatable::xt6
real*8, dimension(: ),allocatable::yt6
real*8, dimension(: ),allocatable::xt7
real*8, dimension(: ),allocatable::yt7
real*8, dimension(: ),allocatable::xt8
real*8, dimension(: ),allocatable::yt8
real*8, dimension(: ),allocatable::xt9
real*8, dimension(: ),allocatable::yt9
real*8, dimension(: ),allocatable::a
real*8, dimension(: ),allocatable::b
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! st
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!***********************************************
!real*8,parameter::light_speed=300000000.0d0
real*8,parameter::jiedian0=8.854187817d-12
real*8,parameter::cidao0=1.256637061d-6
real*8,parameter::pi=3.141592654
integer,parameter::grid_vs_colomn1=10
numbers='0123456789'
light_speed=300000000.0d0
!###############验证程序是用下列
jiedian1=1.
jiedian2=1.
!v_a=1.0d-5
v_a=0.02
!v_b=1.0d-5
v_b=0.02
r=0.25*v_a
write(*,600) '计算开始,此程序是计算波分'
!write(*,600) '当固定波导耦合的长度时,改变频率时,3个不同出口处的透过率'
!open(1,file='c:\ouhe445d.txt',action='write',form='formatted',access='sequential',&
! iostat=ios)
!if (ios/=0) then
!print*,'it is no easy to open the result1.txt,please check again:'
!stop
!end if
do fff=0.415,0.415,0.005
!do fff=0.35,0.45,0.002
!write(1,313) 'fff=',fff
write(*,313) 'fff=',fff
313 format(1x,a,1x,es11.3e2)
!do iii=26,26,1
do iii=50,50,1
length_crystal=v_a*iii
width_crystal= v_b*30
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
goto 8787
write(*,*) ' 矩形2维光子晶体的时域有限差分FORTRAN程序'
write(*,*) ' '
write(*,*) '
write(*,*) '**************************************************'
write(*,111) '请输入你所研究的光子晶体的有关参数:'
write(*,*) ' '
write(*,111) '请输入光子晶体基质的介电常数:(如13.0d0)'
!read(*, *) jiedian1
write(*,111) '请输入光子晶体圆柱空孔的介电常数:(如1.0d0)'
!read(*, *) jiedian2
write(*,111) '请输入光子晶体的长:(如1.0d-4)'
!read(*, *) length_crystal
write(*,111) '请输入光子晶体的宽:(如1.0d-4)'
!read(*, *) width_crystal
write(*,111) '请输入光子晶体X方向的格矢大小:(如1.0d-5)'
!read(*, *) v_a
write(*,111) '请输入光子晶体Y方向的格矢大小: (如1.0d-5)'
!read(*, *) v_b
write(*,111) '请输入光子晶体圆孔的半径:(如 4.0d-6)'
!read(*, *) r
111 format( 1x,A)
!********************************************
8787 nx_colomn=nint(length_crystal/v_a)
nx=nx_colomn*grid_vs_colomn1
dx=v_a/grid_vs_colomn1
ny_colomn=nint(width_crystal/v_b)
if(v_a==v_b) then
grid_vs_colomn2=grid_vs_colomn1
ny=ny_colomn*grid_vs_colomn2
dy=dx
else
grid_vs_colomn2=nint((v_b/v_a)*grid_vs_colomn1)
ny=grid_vs_colomn2*ny_colomn
dy=v_b/grid_vs_colomn2
end if
dt=min(dx,dy)/(2*light_speed)
r_vs_a=r/v_a
r_vs_b=r/v_b
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!*****************************************************************************
zi_jie_shu=0.0
!allocate(jiedian(0:nx,0:ny))
temp1=grid_vs_colomn1/2
temp2=grid_vs_colomn2/2
temp1=5
temp2=5
allocate(jiedian(-temp1:nx+temp1,-temp2:ny+temp2))
allocate(jiedian11(0:nx+2*temp1,0:ny+2*temp2))
do i=0,nx
do j=0,ny
jiedian(i,j)=jiedian1
zi_jie_shu=zi_jie_shu+8
end do
end do
do i=0,nx,grid_vs_colomn1
do j=0,ny,grid_vs_colomn2
jiedian(i,j)=jiedian2
end do
end do
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!***********************************************************************
goto 7878
open(1,file='e:\result2.txt',action='write',form='formatted',access='sequential',&
iostat=ios)
if (ios/=0) then
print*,'it is no easy to open the result1.txt,please check again:'
stop
end if
write(1,*) ' 矩形2维光子晶体的时域有限差分FORTRAN程序'
write(1,*) ' '
write(1,*) '
write(1,*) '**************************************************'
write(1,100) 'physical information of the photonic crystal I code with.&
for example: the width . length of photonic crystal and &
the frequency0 of electromarnetism and lattice constants &
and radius of dielectric colomn &
'
write(1,200) width_crystal
write(1,201) length_crystal
write(1,400) v_a
write(1,401) v_b
write(1,500) r
write(1,501) r_vs_a
write(1,502) r_vs_b
write(*,100) 'physical information of the photonic crystal I code with.&
for example: the width . length of photonic crystal and &
the frequency0 of electromarnetism and lattice constants &
and radius of dielectric colomn &
'
write(*,200) width_crystal
write(*,201) length_crystal
write(*,400) v_a
write(*,401) v_b
write(*,500) r
write(*,501) r_vs_a
write(*,502) r_vs_b
100 format(1x,A)
200 format(1x,' 光子晶体宽:',5x,E13.6E2)
201 format(1x,' 光子晶体长:',5x,E13.6E2)
400 format(1x,' 光子晶体x方向格矢大小:',5x,E13.6E2)
401 format(1x,' 光子晶体y方向格矢大小:',5x,E13.6E2)
500 format(1x,' 光子晶体圆孔半径:',5x,E13.6E2)
501 format(1x,'光子晶体圆孔半径和x方向格矢大小的比值:',5x,E13.6E2)
502 format(1x,'光子晶体圆孔半径和y方向格矢大小的比值:',5x,E13.6E2)
write(1,600) 'messing parameter:'
write(1,700) nx
write(1,701) ny
write(1,800) 'dielectric colomn cells:'
write(1,900) nx_colomn
write(1,901) ny_colomn
write(1,60) dt
write(1,61) dx
write(1,62) dy
write(*,600) 'messing parameter:'
write(*,700) nx
write(*,701) ny
write(*,800) 'dielectric colomn cells:'
write(*,900) nx_colomn
write(*,901) ny_colomn
write(*,60) dt
write(*,61) dx
write(*,62) dy
600 format(1x,A)
700 format(1x,' nx:',5x,I5)
701 format(1x,' ny:',5x,I5)
800 format(1x,A)
900 format(1x,' nx-colomn:',5x,I5)
901 format(1x,' ny-colomn:',5x,I5)
60 format(1x,' timestep dt:',5x,ES13.6E2)
61 format(1x,' spacestep dx:',5x,ES13.6E2)
62 format(1x,' spacestep dy:',5x,ES13.6E2)
write(*,123) '一部分圆孔中心坐标(用来初步检查生成的网格是否正确):'
write(1,123) '一部分圆孔中心坐标(用来初步检查生成的网格是否正确):'
close(1)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
7878 do i=0,3*grid_vs_colomn1
do j=0,3*grid_vs_colomn2
if(abs((jiedian(i,j)-jiedian2))<1.0e-3) then
!write(*,124) i,j
!write(1,124) i,j
else
end if
end do
end do
123 format(1x,A)
124 format(1x,' ','(',i3,',',i3,')')
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!do i=3,6,1
!do j=0,11,1
!write(*,'(1x,ES13.6E2)') even_jiedian(j,3,dx,dy,r,v_a,nx,ny,grid_vs_colomn1,grid_vs_colomn2,jiedian1,jiedian2)
!end do
!end do
!end
!write(*,*) ' 矩形2维光子晶体的10*17网格点处的有效介电常数:'
!write(1,*) ' 矩形2维光子晶体的10*17网格点处的有效介电常数:'
do i=-temp1,nx+temp1,1
do j=-temp2,ny+temp2,1
jiedian(i,j)=even_jiedian(i,j,dx,dy,r,v_a,nx,ny,grid_vs_colomn1,grid_vs_colomn2,jiedian1,jiedian2,nx_colomn)
!if (jiedian(i,j).lt.1.0) write(*,*) jiedian(i,j)
end do
end do
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!quexian
do i=0,nx+2*temp1,1
do j=0,ny+2*temp2,1
temp11=i-temp1
temp22=j-temp2
jiedian11(i,j)=jiedian(temp11,temp22)
end do
end do
nx=nx+2*temp1
ny=ny+2*temp2
! write(*,16) ((jiedian11(i,j),j=0,14),i=72,88)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! write(*,106)
106 format(1x,'经过拓展后的网格参数为:')
!write(*,107) nx
107 format(1x,' x方向的网格单元数nx:',5x,i4)
!write(*,108) ny
108 format(1x,' y方向的网格单元数ny:',5x,i4)
deallocate(jiedian)
!!!!!角点参数的设定
dd=sqrt(real(nx**2+ny**2))
sina1=ny/dd
cosa1=nx/dd
frad=sqrt((dd-1)/dd)
!!!!!角点参数的设定finish
!第二部分
allocate(hx(0:nx,0:ny))
do i=0,nx
do j=0,ny
hx(i,j)=0
zi_jie_shu=zi_jie_shu+8
end do
end do
allocate(hy(0:nx,0:ny))
do i=0,nx
do j=0,ny
hy(i,j)=0
zi_jie_shu=zi_jie_shu+8
end do
end do
allocate(ez(0:nx,0:ny))
do i=0,nx
do j=0,ny
ez(i,j)=0
zi_jie_shu=zi_jie_shu+8
end do
end do
allocate(ez1(0:nx,0:ny))
do i=0,nx
do j=0,ny
ez1(i,j)=0
zi_jie_shu=zi_jie_shu+8
end do
end do
allocate(ez2(0:nx,0:ny))
do i=0,nx
do j=0,ny
ez2(i,j)=0
zi_jie_shu=zi_jie_shu+8
end do
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -