Workbench merge cavity

config/namelist.config

@@ -50,6 +50,8 @@ include_fleapyear=.false.
 
 &run_config
 use_ice=.true.                  ! ocean+ice
+use_cavity=.false.              !
+use_cavity_partial_cell=.false. 
 use_floatice = .false.
 use_sw_pene=.true.
 /

mesh_part/CMakeLists.txt

@@ -19,6 +19,7 @@ add_executable(${PROJECT_NAME} ${sources_Fortran})
 # CMAKE_Fortran_COMPILER_ID will also work if a wrapper is being used (e.g. mpif90 wraps ifort -> compiler id is Intel)
 if(${CMAKE_Fortran_COMPILER_ID} STREQUAL  Intel )
 	target_compile_options(${PROJECT_NAME} PRIVATE -r8 -i4 -fp-model precise -no-prec-div -no-prec-sqrt -fast-transcendentals -xHost -ip)
+#	target_compile_options(${PROJECT_NAME} PRIVATE -r8 -i4 -fp-model precise -no-prec-div -no-prec-sqrt -fast-transcendentals -xHost -ip -g -traceback -check all,noarg_temp_created,bounds,uninit) 
 elseif(${CMAKE_Fortran_COMPILER_ID} STREQUAL  GNU )
 	target_compile_options(${PROJECT_NAME} PRIVATE -fdefault-real-8 -ffree-line-length-none)
 endif()

src/MOD_MESH.F90

@@ -26,8 +26,8 @@ MODULE MOD_MESH
 integer, allocatable, dimension(:,:)        :: elem2D_nodes ! elem2D_nodes(:,n) lists; 3 nodes of element n   
 integer, allocatable, dimension(:,:)        :: edges        ! edge(:,n) lists 2 nodes; edge n
 integer, allocatable, dimension(:,:)        :: edge_tri     ! edge_tri(:,n) lists 2 
-				                           ! elements containing edge n: the first one is to left 
-                                                           ! of the line directed to the second node
+                                                            ! elements containing edge n: the first one is to left 
+                                                            ! of the line directed to the second node
 integer,       allocatable, dimension(:,:)  :: elem_edges   ! elem_edges(:,n) are edges of element n.  
 real(kind=WP), allocatable, dimension(:)    :: elem_area
 real(kind=WP), allocatable, dimension(:,:)  :: edge_dxdy, edge_cross_dxdy
@@ -43,13 +43,47 @@ MODULE MOD_MESH
 real(kind=WP), allocatable, dimension(:,:)  :: gradient_sca ! Coefficients to compute gradient of scalars
                                                            ! on elements
 INTEGER,       ALLOCATABLE, DIMENSION(:)    :: bc_index_nod2D(:)
-                                                           ! vertical structure             
+                                                           ! vertical structure     
+!
+!
+!___vertical mesh info__________________________________________________________                                                           
+! total number of layers
 integer                                     :: nl
+
+! initial layer, mid-depth layer and element depth
 real(kind=WP), allocatable, dimension(:)    :: zbar, Z,elem_depth
-integer,       allocatable, dimension(:)    :: nlevels, nlevels_nod2D
+
+! upper boudnary index of all vertical vertice/element loops, default==1 but when 
+! cavity is used becomes index of cavity-ocean boundary at vertices and elements
+integer,       allocatable, dimension(:)    :: ulevels, ulevels_nod2D, ulevels_nod2D_max
+
+! number of levels at elem and vertices considering bottom topography
+integer,       allocatable, dimension(:)    :: nlevels, nlevels_nod2D, nlevels_nod2D_min
+
+!
+!
+!___horizontal mesh info________________________________________________________
 real(kind=WP), allocatable, dimension(:,:)  :: area, area_inv
 real(kind=WP), allocatable, dimension(:)    :: mesh_resolution
-! Elevation stiffness matrix
+
+!
+!
+!___cavity mesh info____________________________________________________________
+! level index of cavity-ocean boundary at vertices and elements
+! --> see: ulevels, ulevels_nod2D (fvom_main) 
+
+! vertice/element yes=1/no=0 flag if cavity exists
+integer,       allocatable, dimension(:)    :: cavity_flag_n, cavity_flag_e
+
+! depth of cavity-ocean interface
+real(kind=WP), allocatable, dimension(:)    :: cavity_depth
+
+
+real(kind=WP), allocatable, dimension(:,:)  :: cavity_nrst_cavlpnt_xyz
+
+!
+!
+!___Elevation stiffness matrix__________________________________________________
 type(sparse_matrix)                         :: ssh_stiff
 
 !#if defined (__oasis)

src/associate_mesh.h

@@ -21,11 +21,15 @@ real(kind=WP), dimension(:,:), pointer :: gradient_vec
 real(kind=WP), dimension(:,:), pointer :: gradient_sca 
 integer,       dimension(:)  , pointer :: bc_index_nod2D
 real(kind=WP), dimension(:)  , pointer :: zbar, Z, elem_depth
-integer,       dimension(:)  , pointer :: nlevels, nlevels_nod2D
+integer,       dimension(:)  , pointer :: nlevels, nlevels_nod2D, nlevels_nod2D_min
 real(kind=WP), dimension(:,:), pointer :: area, area_inv
 real(kind=WP), dimension(:)  , pointer :: mesh_resolution
 real(kind=WP), dimension(:)  , pointer :: lump2d_north, lump2d_south
 type(sparse_matrix)          , pointer :: ssh_stiff
+!!$integer,       dimension(:)  , pointer :: cavity_lev_nod2D, cavity_lev_elem2D
+integer,       dimension(:)  , pointer :: cavity_flag_n, cavity_flag_e
+real(kind=WP), dimension(:)  , pointer :: cavity_depth
+integer,       dimension(:)  , pointer :: ulevels, ulevels_nod2D, ulevels_nod2D_max
 
 nod2D              => mesh%nod2D              
 elem2D             => mesh%elem2D             
@@ -34,9 +38,9 @@ edge2D_in          => mesh%edge2D_in
 ocean_area         => mesh%ocean_area         
 nl                 => mesh%nl  
 
-!!$coord_nod2D     => mesh%coord_nod2D        
-!!$geo_coord_nod2D => mesh%geo_coord_nod2D    
-!!$elem2D_nodes    => mesh%elem2D_nodes       
+!!$coord_nod2D        => mesh%coord_nod2D        
+!!$geo_coord_nod2D    => mesh%geo_coord_nod2D    
+!!$elem2D_nodes       => mesh%elem2D_nodes       
 !!$edges              => mesh%edges              
 !!$edge_tri           => mesh%edge_tri           
 !!$elem_edges         => mesh%elem_edges         
@@ -59,11 +63,18 @@ nl                 => mesh%nl
 !!$elem_depth         => mesh%elem_depth      
 !!$nlevels            => mesh%nlevels            
 !!$nlevels_nod2D      => mesh%nlevels_nod2D
+!!$nlevels_nod2D_min  => mesh%nlevels_nod2D_min
 !!$area               => mesh%area     
 !!$area_inv           => mesh%area_inv     
 !!$mesh_resolution    => mesh%mesh_resolution    
 !!$ssh_stiff          => mesh%ssh_stiff          
-
+!!$cavity_flag        => mesh%cavity_flag  
+!!$cavity_lev_nod2D   => mesh%cavity_lev_nod2D  
+!!$cavity_lev_elem2D  => mesh%cavity_lev_elem2D  
+!!$cavity_depth       => mesh%cavity_depth  
+!!$ulevels            => mesh%ulevels  
+!!$ulevels_nod2D      => mesh%ulevels_nod2D
+!!$ulevels_nod2D_max  => mesh%ulevels_nod2D_max
 
 coord_nod2D(1:2,1:myDim_nod2D+eDim_nod2D)                  => mesh%coord_nod2D        
 geo_coord_nod2D(1:2,1:myDim_nod2D+eDim_nod2D)              => mesh%geo_coord_nod2D    
@@ -90,9 +101,18 @@ Z(1:mesh%nl-1)                                             => mesh%Z
 elem_depth         => mesh%elem_depth      ! never used, not even allocated
 nlevels(1:myDim_elem2D+eDim_elem2D+eXDim_elem2D)           => mesh%nlevels            
 nlevels_nod2D(1:myDim_nod2D+eDim_nod2D)                    => mesh%nlevels_nod2D
-area(1:mesh%nl,1:myDim_nod2d+eDim_nod2D)                     => mesh%area     
-area_inv(1:mesh%nl,1:myDim_nod2d+eDim_nod2D)                 => mesh%area_inv     
+nlevels_nod2D_min(1:myDim_nod2D+eDim_nod2D)                => mesh%nlevels_nod2D_min
+area(1:mesh%nl,1:myDim_nod2d+eDim_nod2D)                   => mesh%area     
+area_inv(1:mesh%nl,1:myDim_nod2d+eDim_nod2D)               => mesh%area_inv     
 mesh_resolution(1:myDim_nod2d+eDim_nod2D)                  => mesh%mesh_resolution    
 ssh_stiff                                                  => mesh%ssh_stiff
 lump2d_north(1:myDim_nod2d)                                => mesh%lump2d_north
 lump2d_south(1:myDim_nod2d)                                => mesh%lump2d_south
+cavity_flag_n(1:myDim_nod2D+eDim_nod2D)                    => mesh%cavity_flag_n
+cavity_flag_e(1:myDim_elem2D+eDim_elem2D+eXDim_elem2D)     => mesh%cavity_flag_e  
+!!$cavity_lev_nod2D(1:myDim_nod2D+eDim_nod2D)                 => mesh%cavity_lev_nod2D  
+!!$cavity_lev_elem2D(1:myDim_elem2D+eDim_elem2D+eXDim_elem2D) => mesh%cavity_lev_elem2D  
+cavity_depth(1:myDim_nod2D+eDim_nod2D)                     => mesh%cavity_depth  
+ulevels(1:myDim_elem2D+eDim_elem2D+eXDim_elem2D)           => mesh%ulevels            
+ulevels_nod2D(1:myDim_nod2D+eDim_nod2D)                    => mesh%ulevels_nod2D
+ulevels_nod2D_max(1:myDim_nod2D+eDim_nod2D)                => mesh%ulevels_nod2D_max

src/associate_mesh_ini.h

@@ -24,6 +24,8 @@ real(kind=WP), dimension(:)  , pointer :: zbar, Z, elem_depth
 integer,       dimension(:)  , pointer :: nlevels, nlevels_nod2D
 real(kind=WP), dimension(:,:), pointer :: area, area_inv
 real(kind=WP), dimension(:)  , pointer :: mesh_resolution
+integer,       dimension(:)  , pointer :: cavity_flag, ulevels_nod2D, ulevels
+real(kind=WP), dimension(:)  , pointer :: cavity_depth
 type(sparse_matrix)          , pointer :: ssh_stiff
 
 nod2D              => mesh%nod2D              
@@ -62,3 +64,8 @@ area               => mesh%area
 area_inv           => mesh%area_inv     
 mesh_resolution    => mesh%mesh_resolution    
 ssh_stiff          => mesh%ssh_stiff        
+!!$cavity_flag_n      => mesh%cavity_flag_n  
+!!$cavity_flag_e      => mesh%cavity_flag_e
+ulevels_nod2D      => mesh%ulevels_nod2D  
+ulevels            => mesh%ulevels  
+cavity_depth       => mesh%cavity_depth