Skip to content

GitLab

Commit 71b3b061 authored by Knut's avatar Knut
Browse files

save_fluxes: flexible output of velx,vely,velx3d,vely3d

parent ef0aef57
Hide whitespace changes
Inline Side-by-side
Showing with 199 additions and 38 deletions
src/2d/Makefile
View file @ 71b3b061
......@@ -22,6 +22,7 @@ ${LIB}(cfl_check.o) \
${LIB}(bottom_friction.o) \
${LIB}(depth_update.o) \
${LIB}(momentum.o) \
${LIB}(velocity_update.o) \
${LIB}(sealevel.o) \
${LIB}(uv_advect.o) \
${LIB}(uv_diffusion.o) \
......
src/2d/dynamic_allocations_2d.h
View file @ 71b3b061
......@@ -30,6 +30,12 @@
allocate(V(E2DFIELD),stat=rc)
if (rc /= 0) stop 'init_2d: Error allocating memory (V)'
allocate(velx(E2DFIELD),stat=rc)
if (rc /= 0) stop 'init_2d: Error allocating memory (velx)'
allocate(vely(E2DFIELD),stat=rc)
if (rc /= 0) stop 'init_2d: Error allocating memory (vely)'
allocate(UEx(E2DFIELD),stat=rc)
if (rc /= 0) stop 'init_2d: Error allocating memory (UEx)'
......
src/2d/dynamic_declarations_2d.h
View file @ 71b3b061
......@@ -5,6 +5,7 @@
REALTYPE,dimension(:,:),allocatable,target :: t_zo,t_z,D,Dvel,DU,DV
REALTYPE,dimension(:,:),allocatable :: U,V
REALTYPE,dimension(:,:),allocatable :: UEx,VEx
REALTYPE,dimension(:,:),allocatable,target :: velx,vely
REALTYPE,dimension(:,:),allocatable :: ru,rv
REALTYPE,dimension(:,:),allocatable :: Uint,Vint
REALTYPE,dimension(:,:),allocatable :: Uinto,Vinto
......
src/2d/m2d.F90
View file @ 71b3b061
......@@ -43,6 +43,17 @@
REALTYPE,dimension(E2DFIELD),intent(out) :: D,Dvel,DU,DV
end subroutine depth_update
subroutine velocity_update(dt,z,zo,Dvel,U,DU,V,DV,wwm,wwp,missing, &
velx,vely)
use domain, only: imin,imax,jmin,jmax
IMPLICIT NONE
REALTYPE,intent(in) :: dt
REALTYPE,dimension(E2DFIELD),intent(in) :: z,zo,Dvel,U,DU,V,DV
REALTYPE,dimension(E2DFIELD),target,intent(in),optional :: wwm,wwp
REALTYPE,intent(in),optional :: missing
REALTYPE,dimension(E2DFIELD),intent(out) :: velx,vely
end subroutine velocity_update
subroutine uv_advect(U,V,D,Dvel,DU,DV)
use domain, only: imin,imax,jmin,jmax
IMPLICIT NONE
......@@ -412,6 +423,8 @@
end if
call velocity_update(dtm,z,zo,Dvel,U,DU,V,DV,velx=velx,vely=vely)
return
end subroutine postinit_2d
!EOC
......@@ -488,6 +501,7 @@
if (have_boundaries) call update_2d_bdy(loop,bdy2d_ramp)
call sealevel(loop)
call depth_update(zo,z,D,Dvel,DU,DV)
call velocity_update(dtm,z,zo,Dvel,U,DU,V,DV,velx=velx,vely=vely)
if(residual .gt. 0) then
call tic(TIM_INTEGR2D)
......
src/2d/static_2d.h
View file @ 71b3b061
......@@ -8,6 +8,7 @@
REALTYPE,dimension(E2DFIELD),target :: t_zo,t_z,D,Dvel,DU,DV
REALTYPE U(E2DFIELD)
REALTYPE V(E2DFIELD)
REALTYPE,dimension(E2DFIELD),target :: velx,vely
REALTYPE UEx(E2DFIELD)
REALTYPE VEx(E2DFIELD)
REALTYPE ru(E2DFIELD)
......
src/2d/variables_2d.F90
View file @ 71b3b061
......@@ -125,6 +125,7 @@
D = _ZERO_ ; Dvel = _ZERO_
U = _ZERO_; DU = _ZERO_; Uint = _ZERO_; UEx = _ZERO_
V = _ZERO_; DV = _ZERO_; Vint = _ZERO_; VEx = _ZERO_
velx = -9999.0 ; vely = -9999.0
ru = _ZERO_; ruu=_ZERO_; Uinto=_ZERO_
rv = _ZERO_; rvv=_ZERO_; Vinto=_ZERO_
......
src/2d/velocity_update.F90 0 → 100644
View file @ 71b3b061
#include "cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
! !ROUTINE: velocity_update - calculate new velocities.
!
! !INTERFACE:
subroutine velocity_update(dt,z,zo,Dvel,U,DU,V,DV,wwm,wwp,missing, &
velx,vely)
! Note (KK): keep in sync with interface in m2d.F90
!
! !DESCRIPTION:
!
! !USES:
use domain
IMPLICIT NONE
!
! !INPUT PARAMETERS:
REALTYPE,intent(in) :: dt
REALTYPE,dimension(E2DFIELD),intent(in) :: z,zo,Dvel,U,DU,V,DV
REALTYPE,dimension(E2DFIELD),target,intent(in),optional :: wwm,wwp
REALTYPE,intent(in),target,optional :: missing
!
! !OUTPUT PARAMETERS:
REALTYPE,dimension(E2DFIELD),intent(out) :: velx,vely
!
! !OUTPUT PARAMETERS:
!
! !REVISION HISTORY:
! Original author(s): Knut Klingbeil
!
! !LOCAL VARIABLES:
REALTYPE,dimension(E2DFIELD),target :: zeros
REALTYPE,dimension(:,:),pointer :: p_wwm,p_wwp
REALTYPE :: missval
REALTYPE,parameter :: vel_missing=-9999.0
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
integer, save :: Ncall = 0
Ncall = Ncall+1
write(debug,*) 'velocity_update() # ',Ncall
#endif
zeros = _ZERO_
if (present(wwm)) then
p_wwm => wwm
else
p_wwm => zeros
end if
if (present(wwp)) then
p_wwp => wwp
else
p_wwp => zeros
end if
if (present(missing)) then
missval = missing
else
missval = vel_missing
end if
call to_u(imin,jmin,imax,jmax,az, &
dt,grid_type, &
dxv,dyu,arcd1, &
xc,xu,xv,z,zo,Dvel,U,DU,V,DV,p_wwm,p_wwp,missval,velx)
call to_v(imin,jmin,imax,jmax,az, &
dt,grid_type, &
dxv,dyu,arcd1, &
yc,yu,yv,z,zo,Dvel,U,DU,V,DV,p_wwm,p_wwp,missval,vely)
#ifdef DEBUG
write(debug,*) 'Leaving velocity_update()'
write(debug,*)
#endif
return
end subroutine velocity_update
!EOC
!-----------------------------------------------------------------------
! Copyright (C) 2016 - Knut Klingbeil !
!-----------------------------------------------------------------------
src/3d/Makefile
View file @ 71b3b061
......@@ -51,7 +51,8 @@ $(LIB)(hybrid_coordinates.o) \
$(LIB)(bottom_friction_3d.o) \
$(LIB)(uu_momentum_3d.o) \
$(LIB)(vv_momentum_3d.o) \
$(LIB)(ww_momentum_3d.o)
$(LIB)(ww_momentum_3d.o) \
$(LIB)(velocity_update_3d.o)
ifeq ($(GETM_STRUCTURE_FRICTION),true)
OBJ += \
$(LIB)(structure_friction_3d.o)
......
src/3d/dynamic_allocations_3d.h
View file @ 71b3b061
......@@ -7,6 +7,15 @@
allocate(ww(I3DFIELD),stat=rc) ! 3D field for w-velocity
if (rc /= 0) stop 'init_3d: Error allocating memory (ww)'
allocate(velx3d(I3DFIELD),stat=rc) ! 3D field for u-velocity
if (rc /= 0) stop 'init_3d: Error allocating memory (velx3d)'
allocate(vely3d(I3DFIELD),stat=rc) ! 3D field for v-velocity
if (rc /= 0) stop 'init_3d: Error allocating memory (vely3d)'
allocate(w(I3DFIELD),stat=rc) ! 3D field for w-velocity
if (rc /= 0) stop 'init_3d: Error allocating memory (w)'
#ifdef _MOMENTUM_TERMS_
allocate(tdv_u(I3DFIELD),stat=rc) ! 3D field for tdv_u
if (rc /= 0) stop 'init_3d: Error allocating memory (tdv_u)'
......
src/3d/dynamic_declarations_3d.h
View file @ 71b3b061
......@@ -7,6 +7,7 @@
REALTYPE, dimension(:,:,:), allocatable :: uu,vv
REALTYPE, dimension(:,:,:), allocatable, target :: ww
REALTYPE, dimension(:,:,:), allocatable, target :: velx3d,vely3d,w
#ifdef _MOMENTUM_TERMS_
REALTYPE, dimension(:,:,:), allocatable, target :: tdv_u,adv_u,vsd_u
......
src/3d/m3d.F90
View file @ 71b3b061
......@@ -415,6 +415,8 @@
end if
call velocity_update_3d()
return
end subroutine postinit_3d
!EOC
......@@ -632,6 +634,8 @@
end if
call velocity_update_3d()
call slow_terms()
#endif
call tic(TIM_INTEGR3D)
......
src/3d/static_3d.h
View file @ 71b3b061
......@@ -15,6 +15,8 @@
REALTYPE :: uu(I3DFIELD)
REALTYPE :: vv(I3DFIELD)
REALTYPE, target :: ww(I3DFIELD)
REALTYPE,dimension(I3DFIELD),target :: velx3d,vely3d,w
#ifdef _MOMENTUM_TERMS_
REALTYPE :: tdv_u(I3DFIELD)
REALTYPE :: adv_u(I3DFIELD)
......
src/3d/variables_3d.F90
View file @ 71b3b061
......@@ -194,6 +194,7 @@
hn = _ZERO_ ; hvel = _ZERO_ ; hun = _ZERO_ ; hvn = _ZERO_
uu = _ZERO_ ; vv = _ZERO_ ; ww = _ZERO_
velx3d = -9999.0 ; vely3d = -9999.0 ; w = -9999.0
#ifdef _MOMENTUM_TERMS_
tdv_u = _ZERO_ ; adv_u = _ZERO_ ; vsd_u = _ZERO_ ; hsd_u = _ZERO_
......
src/3d/velocity_update_3d.F90 0 → 100644
View file @ 71b3b061
#include "cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
! !ROUTINE: velocity_update_3d - calculate new 3D velocities.
!
! !INTERFACE:
subroutine velocity_update_3d()
!
! !DESCRIPTION:
!
! !USES:
use domain
use m2d, only: velocity_update
use variables_3d, only: kmin
use variables_3d, only: hn,ho,uu,vv,ww,hun,hvn,hvel,velx3d,vely3d,w,dt
IMPLICIT NONE
!
! !INPUT PARAMETERS:
!
! !OUTPUT PARAMETERS:
!
! !REVISION HISTORY:
! Original author(s): Knut Klingbeil
!
! !LOCAL VARIABLES:
integer :: k
REALTYPE,parameter :: vel_missing=-9999.0
!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
integer, save :: Ncall = 0
Ncall = Ncall+1
write(debug,*) 'velocity_update_3d() # ',Ncall
#endif
do k=1,kmax
call velocity_update(dt,hn(:,:,k),ho(:,:,k),hvel(:,:,k), &
uu(:,:,k),hun(:,:,k),vv(:,:,k),hvn(:,:,k), &
wwm=ww(:,:,k-1),wwp=ww(:,:,k), &
velx=velx3d(:,:,k),vely=vely3d(:,:,k))
end do
call to_w(imin,jmin,imax,jmax,kmin,kmax,az, &
dt, &
dxv,dyu,arcd1, &
H,HU,HV,hn,ho,hvel,uu,hun,vv,hvn,ww,vel_missing,w)
#ifdef DEBUG
write(debug,*) 'Leaving velocity_update_3d()'
write(debug,*)
#endif
return
end subroutine velocity_update_3d
!EOC
!-----------------------------------------------------------------------
! Copyright (C) 2016 - Knut Klingbeil !
!-----------------------------------------------------------------------
src/CMakeLists.txt
View file @ 71b3b061
......@@ -184,6 +184,7 @@ add_library(2d OBJECT
2d/uv_diff_2dh.F90
2d/uv_diffusion.F90
2d/variables_2d.F90
2d/velocity_update.F90
)
set_property(TARGET 2d APPEND PROPERTY INCLUDE_DIRECTORIES "${PROJECT_SOURCE_DIR}/2d")
......@@ -231,6 +232,7 @@ add_library(3d OBJECT
3d/uv_advect_3d.F90
3d/uv_diffusion_3d.F90
3d/variables_3d.F90
3d/velocity_update_3d.F90
3d/vv_momentum_3d.F90
3d/ww_momentum_3d.F90
)
......
src/getm/register_all_variables.F90
View file @ 71b3b061
......@@ -315,6 +315,8 @@
call fm%register('D', 'm', 'water depth', standard_name='water depth', fill_value=0.0_rk, data2d=D(_2D_W_), category="2d")
call fm%register('U', 'm2/s', 'transport in local x-direction', standard_name='', data2d=U(_2D_W_), category='2d', output_level=output_level_debug)
call fm%register('V', 'm2/s', 'transport in local y-direction', standard_name='', data2d=V(_2D_W_), category='2d', output_level=output_level_debug)
call fm%register('velx', 'm/s', 'velocity in global x-direction', standard_name='', data2d=velx(_2D_W_), category='2d', fill_value=-9999.0_rk, output_level=output_level_debug)
call fm%register('vely', 'm/s', 'velocity in global y-direction', standard_name='', data2d=vely(_2D_W_), category='2d', fill_value=-9999.0_rk, output_level=output_level_debug)
return
end subroutine register_2d_variables
......@@ -448,6 +450,9 @@
call fm%register('ssvn', 'm', 'elevtion at V-points (3D)', standard_name='', data2d=ssvn(_2D_W_), category='3d', output_level=output_level_debug)
call fm%register('uu', 'm2/s', 'transport in local x-direction (3D)', standard_name='', dimensions=(/id_dim_z/), data3d=uu(_3D_W_), category='3d', output_level=output_level_debug)
call fm%register('vv', 'm2/s', 'transport in local y-direction (3D)', standard_name='', dimensions=(/id_dim_z/), data3d=vv(_3D_W_), category='3d', output_level=output_level_debug)
call fm%register('velx3d', 'm/s', 'velocity in global x-direction (3D)', standard_name='', dimensions=(/id_dim_z/), data3d=velx3d(_3D_W_), category='3d', fill_value=-9999.0_rk, output_level=output_level_debug)
call fm%register('vely3d', 'm/s', 'velocity in global y-direction (3D)', standard_name='', dimensions=(/id_dim_z/), data3d=vely3d(_3D_W_), category='3d', fill_value=-9999.0_rk, output_level=output_level_debug)
call fm%register('w', 'm2/s', 'vertical velocity', standard_name='', dimensions=(/id_dim_z/), data3d=w(_3D_W_), category='3d', fill_value=-9999.0_rk, output_level=output_level_debug)
call fm%register('SS', 's-2', 'shear frequency squared', standard_name='', dimensions=(/id_dim_z/), data3d=SS(_3D_W_), category='3d', output_level=output_level_debug)
end if
......
src/ncdf/save_2d_ncdf.F90
View file @ 71b3b061
......@@ -21,6 +21,7 @@
use domain, only: dxv,dyu,arcd1
use m2d, only: dtm
use variables_2d, only: zo,z,D,Dvel,U,DU,V,DV,res_u,res_v
use variables_2d, only: velx,vely
#ifdef USE_BREAKS
use variables_2d, only: break_stat
#endif
......@@ -93,21 +94,11 @@
! velocites
if (u_id .ne. -1) then
wrk = _ZERO_
call to_u(imin,jmin,imax,jmax,az, &
dtm,grid_type, &
dxv,dyu,arcd1, &
xc,xu,xv,z,zo,Dvel,U,DU,V,DV,wrk,wrk,vel_missing,ws)
err = nf90_put_var(ncid,u_id,ws(_2D_W_),start,edges)
err = nf90_put_var(ncid,u_id,velx(_2D_W_),start,edges)
if (err .NE. NF90_NOERR) go to 10
end if
if (v_id .ne. -1) then
wrk = _ZERO_
call to_v(imin,jmin,imax,jmax,az, &
dtm,grid_type, &
dxv,dyu,arcd1, &
yc,yu,yv,z,zo,Dvel,U,DU,V,DV,wrk,wrk,vel_missing,ws)
err = nf90_put_var(ncid,v_id,ws(_2D_W_),start,edges)
err = nf90_put_var(ncid,v_id,vely(_2D_W_),start,edges)
if (err .NE. NF90_NOERR) go to 10
end if
......
src/ncdf/save_3d_ncdf.F90
View file @ 71b3b061
......@@ -22,6 +22,7 @@
use variables_2d, only: Uinto,Vinto
use variables_3d, only: sseo,ssen,Dn,Dveln,Dun,Dvn
use variables_3d, only: dt,kmin,ho,hn,hvel,uu,hun,vv,hvn,ww,hcc,SS
use variables_3d, only: velx3d,vely3d,w
use variables_3d, only: taubx,tauby
#ifdef _MOMENTUM_TERMS_
use variables_3d, only: tdv_u,adv_u,vsd_u,hsd_u,cor_u,epg_u,ipg_u
......@@ -210,37 +211,15 @@
! velocites
if (uu_id .ne. -1) then
ws(:,:,0) = vel_missing
do k=1,kmax
call to_u(imin,jmin,imax,jmax,az, &
dt,grid_type, &
dxv,dyu,arcd1, &
xc,xu,xv,hn(:,:,k),ho(:,:,k),hvel(:,:,k), &
uu(:,:,k),hun(:,:,k),vv(:,:,k),hvn(:,:,k), &
ww(:,:,k-1),ww(:,:,k),vel_missing,ws(:,:,k))
end do
err = nf90_put_var(ncid,uu_id,ws(_3D_W_),start,edges)
err = nf90_put_var(ncid,uu_id,velx3d(_3D_W_),start,edges)
if (err .NE. NF90_NOERR) go to 10
end if
if (vv_id .ne. -1) then
ws(:,:,0) = vel_missing
do k=1,kmax
call to_v(imin,jmin,imax,jmax,az, &
dt,grid_type, &
dxv,dyu,arcd1, &
yc,yu,yv,hn(:,:,k),ho(:,:,k),hvel(:,:,k), &
uu(:,:,k),hun(:,:,k),vv(:,:,k),hvn(:,:,k), &
ww(:,:,k-1),ww(:,:,k),vel_missing,ws(:,:,k))
end do
err = nf90_put_var(ncid,vv_id,ws(_3D_W_),start,edges)
err = nf90_put_var(ncid,vv_id,vely3d(_3D_W_),start,edges)
if (err .NE. NF90_NOERR) go to 10
end if
if (w_id .ne. -1) then
call to_w(imin,jmin,imax,jmax,kmin,kmax,az, &
dt, &
dxv,dyu,arcd1, &
H,HU,HV,hn,ho,hvel,uu,hun,vv,hvn,ww,vel_missing,ws)
err = nf90_put_var(ncid,w_id,ws(_3D_W_),start,edges)
err = nf90_put_var(ncid,w_id,w(_3D_W_),start,edges)
if (err .NE. NF90_NOERR) go to 10
end if
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment