[Fortran90] how to write increment in dimension

[s_hy]

hi all,

i am new with fortran 90. i have here codes for 2d fdtd polar coordinates. but i have problem handling with increment in dimension in real number (as higlighted below). I have googled how to write the increment in do loops, but do loops only applicable for integer. can anyone help me how to solve this.

thank you

!2d polar coordinates
!reference: Numerical electromagnetic: fdtd method (Umran S Inan, R Marshall)

subroutine test3
implicit none

double precision                   :: f0,miu0,eps0
double precision                   :: delta, dt,dr,dp
double precision                   :: vp,lam0,rph,rmh           !lam0 = lambda          
integer                                 :: pp
integer                                 :: rr,m,n
integer                                 :: t
integer                                 :: i,j
double precision,dimension(202,202) :: x,y
double precision,dimension(202,202) :: Ez,Hp,Hr    !Ez,Hy,Hx
double precision                   :: Ein 
double precision, parameter  :: pi = 3.14159265
double precision                   :: r
double precision                   :: p  
character(len=20)                 :: filename


 f0   = 30.e9
 miu0 = 4.0*pi*1.0e-7
 eps0 = 8.854e-12
 vp   = 1/sqrt(eps0*miu0)
 lam0 = vp/f0
 dr   = 0.05*lam0
 dt   = dr/vp
 
! space grid 201x201      
 rr = 201
[COLOR="Blue"] r = ,dr*(rr-1),dr[/COLOR]

!define deltaphi so that there are 200 points in 2pi rad
 dp = 2*pi/200        !in rad

[COLOR="Blue"] p = 1,2*pi-dp,dp[/COLOR]
pp = size(p)

!initialize Ez,Hp,Hr to zero at t=0
do m = rr,pp
 do n = rr,pp
  Hp(m,n) = 0
  Hr(m,n) = 0
  Ez(m,n) = 0
 end do
end do

! coodinate transformation
do m = rr,pp
 do n = rr,pp
  x(m,n) = 0.0
  y(m,n) = 0.0
 end do
end do

do m = 1,rr
 do n = 1,pp
  x(m,n) = r(m)*cos(p(n))
  y(m,n) = r(m)*sin(p(n))
 end do
end do



!timesteps
do t = 1,200
   
  write (filename, "('data',I3.3,'.dat')") n
  open (unit=130,file=filename)
 
!initiate sinusoidal wavepulse at center
!  Ein  = 0.0
!  Ein(t) = sin(2*pi*f0*t*dt)
  Ez(1,1) = sin(2*pi*f0*t*dt)
  !print *, 'Ein(t)=',Ein(t)

!update magnetic field equation for Hr
  do m = 2,rr
   do n = 1,pp-1
    Hr(m,n) = Hr(m,n) - ((dt/miu0/r(m)/dp)*(Ez(m,n+1) - Ez(m,n)))
    !print *,'i=',i,'j=',j,'Hr(i,j)=',Hr(i,j)
   end do
  end do

! correction for 2pi continuity
  do m = 2,rr
   Hr(m,pp) = Hr(m,pp) - ((dt/miu0/r(m)/dp)*(Ez(m,1) - Ez(m,pp)))
  end do

!update magnetic field for Hp
  do m = 1,rr-1
   do n = 1,pp
    Hp(m,n) = Hp(m,n) + ((dt/miu0/dr)*(Ez(m+1,n) - Ez(m,n)))
    !print *,'m=',m,'n=',n,'Hp(m,n)=',Hp(m,n)
   end do
  end do

!fix end point where phi = 2pi
  do m = 2,rr
    Ez(m,1) = Ez(m,1)+(dt/eps0/r(m)/dr)*(rph*Hp(m,1)-rmh*Hp(m-1,1))&
-((dt/eps0/r(m)/dp)*(Hr(m,1)-Hr(m,pp)))
  end do

!update electric field equation 
  do m = 2,rr
   do n = 2,pp
    rph = r(m)+dr/2
    rmh = r(m)-dr/2   
    Ez(m,n) = Ez(m,n)+(dt/eps0/r(m)/dr)*(rph*Hp(m,n)-rmh*Hp(m-1,n))&
-(dt/eps0/r(m)/dp*(Hr(m,n)-Hr(m,n-1)))
!    write (130,*) m,n,Ez(m,n)
!    if (m == 200) write (130,*) ' '
   end do
  end do

close (unit=130) 

end do !n

end SUBROUTINE test3

 

[SteamKing]

Your assignment statements in blue for r = and p = don't make sense in Fortran.

If you want to create a series of real values for r and p which are equally spaced,
then define r and p as real arrays. You already know the size of the array you will need given the
number of points you want to create in the interval. Arrays in F90 can be dimensioned in a flexible fashion
depending on your requirements.

See: http://www.cs.mtu.edu/~shene/COURSES/cs201/NOTES/F90-Array.pdf

You can also use the implied do loop to initialize the values of your arrays.

 

[s_hy]

actually i have MATLAB code and trying to write in fortran. however i have problem with the interval part for r and p. however, for r i have write it id do loops

!increment r = 1,dr*(rr-1),dr
do ir = 1,rr-1
  r(ir) = (ir-1)*dr
  !print *, 'ir=',ir,'r=',r
end do

and it is work. but i have problem with p, which is for p the last value should be 2*pi-dp, with the intervalor increment is dp...and size for p is 200. ( p=0:2*pi-dp:dp). can anyone teach me how to write this increment using do loop.

thank you

 

[gsal]

n  = 200
dp = 2*pi/n
do ip = 0,n-1
   p = ip*dp
   .
   .
   .
end do

By the way, your loop in the previous post seems short one extra dr...by having ir go from 1 to (rr-1) and additionally multiplying dr times (ir-1)...you have one extra -1; also, you may not want to go from 1, but 0 as my example shows to go from 0 to 2*pi-dp

 

[s_hy]

thank you gsal. it's work perfectly.