Commit 1a6cf3c3 authored by Bjarne Buchmann's avatar Bjarne Buchmann
Browse files

Update of some constants from float

parent 8eba85d9
...@@ -102,15 +102,15 @@ ...@@ -102,15 +102,15 @@
if ((((i.eq.1).or.(i.eq.imax)).and.(j.ge.1).and.(j.le.jmax)).or. & if ((((i.eq.1).or.(i.eq.imax)).and.(j.ge.1).and.(j.le.jmax)).or. &
(((j.eq.1).or.(j.eq.jmax)).and.(i.ge.1).and.(i.le.imax))) & (((j.eq.1).or.(j.eq.jmax)).and.(i.ge.1).and.(i.le.imax))) &
z(i,j)=(1.-kk)*z(i,j) z(i,j)=(1.-kk)*z(i,j)
kk=0.5625 kk=0.5625d0
if ((((i.eq.2).or.(i.eq.imax-1)).and.(j.ge.2).and.(j.le.jmax-1)).or. & if ((((i.eq.2).or.(i.eq.imax-1)).and.(j.ge.2).and.(j.le.jmax-1)).or. &
(((j.eq.2).or.(j.eq.jmax-1)).and.(i.ge.2).and.(i.le.imax-1))) & (((j.eq.2).or.(j.eq.jmax-1)).and.(i.ge.2).and.(i.le.imax-1))) &
z(i,j)=(1.-kk)*z(i,j) z(i,j)=(1.-kk)*z(i,j)
kk=0.25 kk=0.25d0
if ((((i.eq.3).or.(i.eq.imax-2)).and.(j.ge.3).and.(j.le.jmax-2)).or. & if ((((i.eq.3).or.(i.eq.imax-2)).and.(j.ge.3).and.(j.le.jmax-2)).or. &
(((j.eq.3).or.(j.eq.jmax-2)).and.(i.ge.3).and.(i.le.imax-2))) & (((j.eq.3).or.(j.eq.jmax-2)).and.(i.ge.3).and.(i.le.imax-2))) &
z(i,j)=(1.-kk)*z(i,j) z(i,j)=(1.-kk)*z(i,j)
kk=0.0625 kk=0.0625d0
if ((((i.eq.4).or.(i.eq.imax-3)).and.(j.ge.4).and.(j.le.jmax-3)).or. & if ((((i.eq.4).or.(i.eq.imax-3)).and.(j.ge.4).and.(j.le.jmax-3)).or. &
(((j.eq.4).or.(j.eq.jmax-3)).and.(i.ge.4).and.(i.le.imax-3))) & (((j.eq.4).or.(j.eq.jmax-3)).and.(i.ge.4).and.(i.le.imax-3))) &
z(i,j)=(1.-kk)*z(i,j) z(i,j)=(1.-kk)*z(i,j)
...@@ -118,7 +118,7 @@ ...@@ -118,7 +118,7 @@
#endif #endif
#ifdef USE_BREAKS #ifdef USE_BREAKS
if (z(i,j)+H(i,j) .lt. 0.9*min_depth .and. & if (z(i,j)+H(i,j) .lt. 0.9d0*min_depth .and. &
break_mask(i,j) .eq. 0 ) then break_mask(i,j) .eq. 0 ) then
break_mask(i,j)=1 break_mask(i,j)=1
break_stat(i,j)=break_stat(i,j)+1 break_stat(i,j)=break_stat(i,j)+1
......
...@@ -47,8 +47,8 @@ ...@@ -47,8 +47,8 @@
#endif #endif
integer, public, parameter :: UPSTREAM=1,UPSTREAM_SPLIT=2,P2=3 integer, public, parameter :: UPSTREAM=1,UPSTREAM_SPLIT=2,P2=3
integer, public, parameter :: Superbee=4,MUSCL=5,P2_PDM=6,FCT=7 integer, public, parameter :: Superbee=4,MUSCL=5,P2_PDM=6,FCT=7
REALTYPE, public, parameter :: one6th=1./6. REALTYPE, public, parameter :: one6th=_ONE_/6
REALTYPE, public, parameter :: ONE=_ONE_,TWO=2.*_ONE_ REALTYPE, public, parameter :: ONE=_ONE_,TWO=_TWO_
! !
! !REVISION HISTORY: ! !REVISION HISTORY:
! Original author(s): Karsten Bolding & Hans Burchard ! Original author(s): Karsten Bolding & Hans Burchard
......
...@@ -56,7 +56,7 @@ ...@@ -56,7 +56,7 @@
character(len=PATH_MAX) :: bdyfile_3d character(len=PATH_MAX) :: bdyfile_3d
REALTYPE :: ip_fac=_ONE_ REALTYPE :: ip_fac=_ONE_
integer :: vel_check=0 integer :: vel_check=0
REALTYPE :: min_vel=-4.,max_vel=4. REALTYPE :: min_vel=-4*_ONE_,max_vel=4*_ONE_
! !
! !REVISION HISTORY: ! !REVISION HISTORY:
! Original author(s): Karsten Bolding & Hans Burchard ! Original author(s): Karsten Bolding & Hans Burchard
...@@ -249,8 +249,8 @@ ...@@ -249,8 +249,8 @@
num=avmback num=avmback
nuh=avhback nuh=avhback
#else #else
num=1.e-15 num=1.d-15
nuh=1.e-15 nuh=1.d-15
#endif #endif
! Needed for interpolation of temperature and salinity ! Needed for interpolation of temperature and salinity
......
...@@ -33,12 +33,12 @@ ...@@ -33,12 +33,12 @@
character(len=PATH_MAX) :: salt_file="t_and_s.nc" character(len=PATH_MAX) :: salt_file="t_and_s.nc"
integer :: salt_field_no=1 integer :: salt_field_no=1
character(len=32) :: salt_name='salt' character(len=32) :: salt_name='salt'
REALTYPE :: salt_const=35. REALTYPE :: salt_const=35*_ONE_
integer :: salt_hor_adv=1,salt_ver_adv=1 integer :: salt_hor_adv=1,salt_ver_adv=1
integer :: salt_adv_split=0 integer :: salt_adv_split=0
REALTYPE :: salt_AH=-_ONE_ REALTYPE :: salt_AH=-_ONE_
integer :: salt_check=0 integer :: salt_check=0
REALTYPE :: min_salt=0.,max_salt=40. REALTYPE :: min_salt=_ZERO_,max_salt=40*_ONE_
! !
! !REVISION HISTORY: ! !REVISION HISTORY:
! Original author(s): Karsten Bolding & Hans Burchard ! Original author(s): Karsten Bolding & Hans Burchard
...@@ -265,7 +265,8 @@ salt_field_no=1 ...@@ -265,7 +265,8 @@ salt_field_no=1
S = _ZERO_ S = _ZERO_
do i=1,100 do i=1,100
do k=1,kmax do k=1,kmax
S(i,2,k)=30.*(1.- tanh(float(i-1)*0.05)) ! S(i,2,k)=30.*(1.- tanh(float(i-1)*0.05))
S(i,2,k)=(30*_ONE_)*(_ONE_- tanh((i-1)*_ONE_/20))
end do end do
end do end do
#endif #endif
...@@ -277,22 +278,22 @@ salt_field_no=1 ...@@ -277,22 +278,22 @@ salt_field_no=1
#ifdef ARKONA_TEST #ifdef ARKONA_TEST
do i=100,135 do i=100,135
do j=256,257 do j=256,257
if (az(i,j).ge.1) S(i,j,0:kmax) = 25. if (az(i,j).ge.1) S(i,j,0:kmax) = 25*_ONE_
end do end do
end do end do
do i=26,27 do i=26,27
do j=77,100 do j=77,100
S(i,j,0:kmax) = 8. S(i,j,0:kmax) = 8*_ONE_
end do end do
end do end do
#endif #endif
#ifdef INTERLEAVING_TEST #ifdef INTERLEAVING_TEST
S(2:6,2,1:20) = 12. S(2:6,2,1:20) = 12*_ONE_
#endif #endif
#ifdef SLOPE_TEST #ifdef SLOPE_TEST
do i=81,82 do i=81,82
do j=42,43 do j=42,43
S(i,j,0:kmax) = 25. S(i,j,0:kmax) = 25*_ONE_
end do end do
end do end do
#endif #endif
...@@ -300,17 +301,17 @@ salt_field_no=1 ...@@ -300,17 +301,17 @@ salt_field_no=1
j=2 j=2
if (imax.eq.102) then if (imax.eq.102) then
do i=2,21 do i=2,21
S(i,j,0:kmax) = 25. S(i,j,0:kmax) = 25*_ONE_
end do end do
end if end if
if (imax.eq.302) then if (imax.eq.302) then
do i=2,61 do i=2,61
S(i,j,0:kmax) = 25. S(i,j,0:kmax) = 25*_ONE_
end do end do
end if end if
if (imax.eq.902) then if (imax.eq.902) then
do i=2,181 do i=2,181
S(i,j,0:kmax) = 25. S(i,j,0:kmax) = 25*_ONE_
end do end do
end if end if
#endif #endif
...@@ -469,11 +470,12 @@ salt_field_no=1 ...@@ -469,11 +470,12 @@ salt_field_no=1
j=2 j=2
do k=1,kmax do k=1,kmax
do i=1,100 do i=1,100
kk= 1.- tanh(float(i-1)*0.05) ! kk= 1.- tanh(float(i-1)*0.05)
S(i,j,k)=(1.-kk)*S(i,j,k)+kk*SRelax kk= _ONE_- tanh((i-1)*_ONE_/20)
S(i,j,k)=(_ONE_-kk)*S(i,j,k)+kk*SRelax
end do end do
end do end do
S(imax-1,2,:)=0. !river S(imax-1,2,:)=_ZERO_ !river
#endif #endif
...@@ -512,9 +514,9 @@ salt_field_no=1 ...@@ -512,9 +514,9 @@ salt_field_no=1
if (kmax.gt.1) then if (kmax.gt.1) then
! Auxilury terms, old and new time level, ! Auxilury terms, old and new time level,
do k=1,kmax-1 do k=1,kmax-1
auxo(k)=2*(1-cnpar)*dt*(nuh(i,j,k)+avmols)/ & auxo(k)=_TWO_*(1-cnpar)*dt*(nuh(i,j,k)+avmols)/ &
(hn(i,j,k+1)+hn(i,j,k)) (hn(i,j,k+1)+hn(i,j,k))
auxn(k)=2* cnpar *dt*(nuh(i,j,k)+avmols)/ & auxn(k)=_TWO_* cnpar *dt*(nuh(i,j,k)+avmols)/ &
(hn(i,j,k+1)+hn(i,j,k)) (hn(i,j,k+1)+hn(i,j,k))
end do end do
...@@ -574,12 +576,12 @@ salt_field_no=1 ...@@ -574,12 +576,12 @@ salt_field_no=1
#ifdef ARKONA_TEST #ifdef ARKONA_TEST
do i=100,135 do i=100,135
do j=256,257 do j=256,257
if (az(i,j).ge.1) S(i,j,0:kmax) = 25. if (az(i,j).ge.1) S(i,j,0:kmax) = 25*_ONE_
end do end do
end do end do
do i=26,27 do i=26,27
do j=77,100 do j=77,100
S(i,j,0:kmax) = 8. S(i,j,0:kmax) = 8*_ONE_
end do end do
end do end do
#endif #endif
...@@ -587,7 +589,7 @@ salt_field_no=1 ...@@ -587,7 +589,7 @@ salt_field_no=1
#ifdef SLOPE_TEST #ifdef SLOPE_TEST
do i=81,82 do i=81,82
do j=42,43 do j=42,43
S(i,j,0:kmax) = 25. S(i,j,0:kmax) = 25*_ONE_
end do end do
end do end do
#endif #endif
......
...@@ -220,10 +220,10 @@ ...@@ -220,10 +220,10 @@
end if end if
fc=f(i ,j,k) ! central fc=f(i ,j,k) ! central
fd=f(i+1,j,k) ! downstream fd=f(i+1,j,k) ! downstream
if (abs(fd-fc) .gt. 1e-10) then if (abs(fd-fc) .gt. 1.d-10) then
r=(fc-fu)/(fd-fc) ! slope ratio r=(fc-fu)/(fd-fc) ! slope ratio
else else
r=(fc-fu)*1.e10 r=(fc-fu)*1.d10
end if end if
else else
c=-uu(i,j,k)/hun(i,j,k)*dt/delxu(i,j) c=-uu(i,j,k)/hun(i,j,k)*dt/delxu(i,j)
...@@ -234,10 +234,10 @@ ...@@ -234,10 +234,10 @@
end if end if
fc=f(i+1,j,k) ! central fc=f(i+1,j,k) ! central
fd=f(i ,j,k) ! downstream fd=f(i ,j,k) ! downstream
if (abs(fc-fd) .gt. 1e-10) then if (abs(fc-fd) .gt. 1.d-10) then
r=(fu-fc)/(fc-fd) ! slope ratio r=(fu-fc)/(fc-fd) ! slope ratio
else else
r=(fu-fc)*1.e10 r=(fu-fc)*1.d10
end if end if
end if end if
select case (method) select case (method)
...@@ -248,12 +248,12 @@ ...@@ -248,12 +248,12 @@
limit=Phi limit=Phi
else else
limit=max(_ZERO_,min(Phi,_TWO_/(_ONE_-c), & limit=max(_ZERO_,min(Phi,_TWO_/(_ONE_-c), &
_TWO_*r/(c+1.e-10))) _TWO_*r/(c+1.d-10)))
end if end if
case (Superbee) case (Superbee)
limit=max(_ZERO_,min(_ONE_,_TWO_*r),min(r,_TWO_)) limit=max(_ZERO_,min(_ONE_,_TWO_*r),min(r,_TWO_))
case (MUSCL) case (MUSCL)
limit=max(_ZERO_,min(_TWO_,_TWO_*r,0.5*(_ONE_+r))) limit=max(_ZERO_,min(_TWO_,_TWO_*r,_HALF_*(_ONE_+r)))
case default case default
FATAL 'Not so good - do_advection_3d()' FATAL 'Not so good - do_advection_3d()'
stop 'u_split_adv' stop 'u_split_adv'
......
...@@ -145,7 +145,7 @@ ...@@ -145,7 +145,7 @@
else else
cu(i,j,k)=uu(i,j,k)*f(i+1,j,k) cu(i,j,k)=uu(i,j,k)*f(i+1,j,k)
end if end if
if ((AH.gt.0.).and.(az(i,j).gt.0).and.(az(i+1,j).gt.0)) & if ((AH.gt._ZERO_).and.(az(i,j).gt.0).and.(az(i+1,j).gt.0)) &
cu(i,j,k)=cu(i,j,k)-AH*(f(i+1,j,k)-f(i,j,k))/delxu(i,j) & cu(i,j,k)=cu(i,j,k)-AH*(f(i+1,j,k)-f(i,j,k))/delxu(i,j) &
*_HALF_*(hn(i+1,j,k)+hn(i,j,k)) *_HALF_*(hn(i+1,j,k)+hn(i,j,k))
end do end do
...@@ -176,7 +176,7 @@ ...@@ -176,7 +176,7 @@
else else
cu(i,j,k)=vv(i,j,k)*f(i,j+1,k) cu(i,j,k)=vv(i,j,k)*f(i,j+1,k)
end if end if
if ((AH.gt.0.).and.(az(i,j).gt.0).and.(az(i,j+1).gt.0)) & if ((AH.gt._ZERO_).and.(az(i,j).gt.0).and.(az(i,j+1).gt.0)) &
cu(i,j,k)=cu(i,j,k)-AH*(f(i,j+1,k)-f(i,j,k))/delyv(i,j) & cu(i,j,k)=cu(i,j,k)-AH*(f(i,j+1,k)-f(i,j,k))/delyv(i,j) &
*_HALF_*(hn(i,j+1,k)+hn(i,j,k)) *_HALF_*(hn(i,j+1,k)+hn(i,j,k))
end do end do
......
...@@ -137,7 +137,7 @@ ...@@ -137,7 +137,7 @@
#ifdef DEBUG #ifdef DEBUG
integer, save :: Ncall = 0 integer, save :: Ncall = 0
Ncall = Ncall+1 Ncall = Ncall+1
write(debug,*) 'D3uvDiff() # ',Ncall write(debug,*) 'uv_diffusion_3d() # ',Ncall
#endif #endif
call tic(TIM_UVDIFF3D) call tic(TIM_UVDIFF3D)
...@@ -152,7 +152,7 @@ ...@@ -152,7 +152,7 @@
PP(i,j,k)=_ZERO_ PP(i,j,k)=_ZERO_
if (az(i,j) .ge. 1) then if (az(i,j) .ge. 1) then
if (k .ge. kumin(i,j)) then if (k .ge. kumin(i,j)) then
PP(i,j,k)=2.*Am*DYC*hn(i,j,k) & PP(i,j,k)=_TWO_*Am*DYC*hn(i,j,k) &
*(uu(i,j,k)/hun(i,j,k)-uu(i-1,j,k)/hun(i-1,j,k))/DXC *(uu(i,j,k)/hun(i,j,k)-uu(i-1,j,k)/hun(i-1,j,k))/DXC
end if end if
end if end if
...@@ -219,7 +219,7 @@ ...@@ -219,7 +219,7 @@
PP(i,j,k)=_ZERO_ PP(i,j,k)=_ZERO_
if (ax(i,j) .ge. 1) then if (ax(i,j) .ge. 1) then
if (k .ge. kumin(i,j)) then if (k .ge. kumin(i,j)) then
PP(i,j,k)=Am*0.5*(hvn(i+1,j,k)+hvn(i,j,k))*DXX & PP(i,j,k)=Am*_HALF_*(hvn(i+1,j,k)+hvn(i,j,k))*DXX &
*((uu(i,j+1,k)/hun(i,j+1,k)-uu(i,j,k)/hun(i,j,k))/DYX & *((uu(i,j+1,k)/hun(i,j+1,k)-uu(i,j,k)/hun(i,j,k))/DYX &
+(vv(i+1,j,k)/hvn(i+1,j,k)-vv(i,j,k)/hvn(i,j,k))/DXX ) +(vv(i+1,j,k)/hvn(i+1,j,k)-vv(i,j,k)/hvn(i,j,k))/DXX )
end if end if
...@@ -252,7 +252,7 @@ ...@@ -252,7 +252,7 @@
do i=imin,imax ! PP defined on T-points do i=imin,imax ! PP defined on T-points
if (az(i,j) .ge. 1) then if (az(i,j) .ge. 1) then
if (k .ge. kvmin(i,j)) then if (k .ge. kvmin(i,j)) then
PP(i,j,k)=2.*Am*DXC*hn(i,j,k) & PP(i,j,k)=_TWO_*Am*DXC*hn(i,j,k) &
*(vv(i,j,k)/hvn(i,j,k)-vv(i,j-1,k)/hvn(i,j-1,k))/DYC *(vv(i,j,k)/hvn(i,j,k)-vv(i,j-1,k)/hvn(i,j-1,k))/DYC
end if end if
end if end if
......
...@@ -128,10 +128,10 @@ ...@@ -128,10 +128,10 @@
end if end if
fc=f(i,j ,k) ! central fc=f(i,j ,k) ! central
fd=f(i,j+1,k) ! downstream fd=f(i,j+1,k) ! downstream
if (abs(fd-fc) .gt. 1e-10) then if (abs(fd-fc) .gt. 1.d-10) then
r=(fc-fu)/(fd-fc) ! slope ratio r=(fc-fu)/(fd-fc) ! slope ratio
else else
r=(fc-fu)*1.e10 r=(fc-fu)*1.d10
end if end if
else else
c=-vv(i,j,k)/hvn(i,j,k)*dt/delyv(i,j) c=-vv(i,j,k)/hvn(i,j,k)*dt/delyv(i,j)
...@@ -142,10 +142,10 @@ ...@@ -142,10 +142,10 @@
end if end if
fc=f(i,j+1,k) ! central fc=f(i,j+1,k) ! central
fd=f(i,j ,k) ! downstream fd=f(i,j ,k) ! downstream
if (abs(fc-fd) .gt. 1e-10) then if (abs(fc-fd) .gt. 1.d-10) then
r=(fu-fc)/(fc-fd) ! slope ratio r=(fu-fc)/(fc-fd) ! slope ratio
else else
r=(fu-fc)*1.e10 r=(fu-fc)*1.d10
end if end if
end if end if
select case (method) select case (method)
...@@ -155,7 +155,7 @@ ...@@ -155,7 +155,7 @@
if (method.eq.P2) then if (method.eq.P2) then
limit=Phi limit=Phi
else else
limit=max(_ZERO_,min(Phi,_TWO_/(_ONE_-c),_TWO_*r/(c+1.e-10))) limit=max(_ZERO_,min(Phi,_TWO_/(_ONE_-c),_TWO_*r/(c+1.d-10)))
end if end if
case (Superbee) case (Superbee)
limit=max(_ZERO_, min(_ONE_,_TWO_*r), min(r,_TWO_) ) limit=max(_ZERO_, min(_ONE_,_TWO_*r), min(r,_TWO_) )
......
...@@ -117,10 +117,10 @@ ...@@ -117,10 +117,10 @@
end if end if
fc=f(i,j,k ) ! central fc=f(i,j,k ) ! central
fd=f(i,j,k+1) ! downstream fd=f(i,j,k+1) ! downstream
if (abs(fd-fc) .gt. 1e-10) then if (abs(fd-fc) .gt. 1.d-10) then
r=(fc-fu)/(fd-fc) ! slope ratio r=(fc-fu)/(fd-fc) ! slope ratio
else else
r=(fc-fu)*1.e10 r=(fc-fu)*1.d10
end if end if
else else
c=-ww(i,j,k)*dt/(_HALF_*(hi(i,j,k)+hi(i,j,k+1))) c=-ww(i,j,k)*dt/(_HALF_*(hi(i,j,k)+hi(i,j,k+1)))
...@@ -131,10 +131,10 @@ ...@@ -131,10 +131,10 @@
end if end if
fc=f(i,j,k+1) ! central fc=f(i,j,k+1) ! central
fd=f(i,j,k ) ! downstream fd=f(i,j,k ) ! downstream
if (abs(fc-fd) .gt. 1e-10) then if (abs(fc-fd) .gt. 1.d-10) then
r=(fu-fc)/(fc-fd) ! slope ratio r=(fu-fc)/(fc-fd) ! slope ratio
else else
r=(fu-fc)*1.e10 r=(fu-fc)*1.d10
end if end if
end if end if
select case (method) select case (method)
...@@ -144,7 +144,7 @@ ...@@ -144,7 +144,7 @@
if (method.eq.P2) then if (method.eq.P2) then
limit=Phi limit=Phi
else else
limit=max(_ZERO_,min(Phi,_TWO_/(_ONE_-c),_TWO_*r/(c+1.e-10))) limit=max(_ZERO_,min(Phi,_TWO_/(_ONE_-c),_TWO_*r/(c+1.d-10)))
end if end if
case (Superbee) case (Superbee)
limit=max(_ZERO_, min(_ONE_, _TWO_*r), min(r,_TWO_) ) limit=max(_ZERO_, min(_ONE_, _TWO_*r), min(r,_TWO_) )
......
...@@ -170,10 +170,10 @@ ...@@ -170,10 +170,10 @@
else else
fd=f(i,j,k) ! downstream fd=f(i,j,k) ! downstream
end if end if
if (abs(fd-fc) .gt. 1e-10) then if (abs(fd-fc) .gt. 1.d-10) then
r=(fc-fu)/(fd-fc) ! slope ratio r=(fc-fu)/(fd-fc) ! slope ratio
else else
r=(fc-fu)*1.e10 r=(fc-fu)*1.d10
end if end if
else else
if (k.lt.kmax) then if (k.lt.kmax) then
...@@ -197,10 +197,10 @@ ...@@ -197,10 +197,10 @@
fc=f(i,j,k) ! central fc=f(i,j,k) ! central
end if end if
fd=f(i,j,k ) ! downstream fd=f(i,j,k ) ! downstream
if (abs(fc-fd) .gt. 1e-10) then if (abs(fc-fd) .gt. 1.d-10) then
r=(fu-fc)/(fc-fd) ! slope ratio r=(fu-fc)/(fc-fd) ! slope ratio
else else
r= (fu-fc)*1.e10 r= (fu-fc)*1.d10
end if end if
end if end if
select case (method) select case (method)
...@@ -210,7 +210,7 @@ ...@@ -210,7 +210,7 @@
if (method.eq.P2) then