From 36275aed478d21de3c9051cf7155a2a0a50ea2b7 Mon Sep 17 00:00:00 2001 From: Alistair Adcroft Date: Tue, 13 Sep 2016 14:47:36 -0400 Subject: [PATCH] *Altered initial conditions for circle_obcs - The circle_obcs experiment ostensibly tests the Flather open boundary conditions except it was modified to use 10 layers and also test the Orlanski boundary conditions too. The initialization code appears to have been written for a 2-layer system and essentially only forces a very high-mode. This change produces a mixed mode initial condition (thickness anomaly linear with z) but with most energy in the first baroclinic mode. - Changes circle_obcs experiment only. --- src/user/circle_obcs_initialization.F90 | 12 ++++++++++-- 1 file changed, 10 insertions(+), 2 deletions(-) diff --git a/src/user/circle_obcs_initialization.F90 b/src/user/circle_obcs_initialization.F90 index 80e326674f..e16cb73ae5 100644 --- a/src/user/circle_obcs_initialization.F90 +++ b/src/user/circle_obcs_initialization.F90 @@ -93,8 +93,16 @@ subroutine circle_obcs_initialize_thickness(h, G, GV, param_file) ! if (rad <= 6.*diskrad) h(i,j,k) = h(i,j,k)+10.0*exp( -0.5*( rad**2 ) ) rad = min( rad, 1. ) ! Flatten outside radius of diskrad rad = rad*(2.*asin(1.)) ! Map 0-1 to 0-pi - h(i,j,k) = h(i,j,k) + 10.0*0.5*(1.+cos(rad)) ! cosine bell - h(i,j,k-1) = h(i,j,k-1) - 10.0*0.5*(1.+cos(rad)) ! cosine bell + if (Nz==1) then + ! The model is barotropic + h(i,j,k) = h(i,j,k) + 1.0*0.5*(1.+cos(rad)) ! cosine bell + else + ! The model is baroclinic + do k = 1, Nz + h(i,j,k) = h(i,j,k) - 0.5*(1.+cos(rad)) & ! cosine bell + * 5.0 * real( 2*k-Nz ) + enddo + endif enddo ; enddo end subroutine circle_obcs_initialize_thickness