Wombat

.travis.yml

@@ -19,6 +19,8 @@ env:
         - TYPE=EBM
         - TYPE=ACCESS-OM
         - TYPE=ACCESS-CM
+        - TYPE=ACCESS-OM-BGC
+        - TYPE=ACCESS-ESM
     global:
         - FC=gfortran-4.8
         - OMPI_FC=${FC}

exp/MOM_compile.csh

@@ -54,8 +54,10 @@ if ( $help ) then
     echo "             EBM       : ocean-seaice-land-atmosphere coupled model with energy balance atmosphere"
     echo "             ACCESS-CM : ocean component of ACCESS-CM model."
     echo "             ACCESS-OM : ocean component of ACCESS-OM model."
+    echo "             ACCESS-ESM : ocean component of ACCESS-ESM model with CSIRO BGC  (Wombat)."
+    echo "             ACCESS-OM-BGC: ocean component of ACCESS-OM model with CSIRO BGC  (Wombat)."
     echo
-    echo "--platform   followed by the platform name that has a corresponding environs file in the ../bin dir, default is gfortran"
+    echo "--platform   followed by the platform name that has a corresponding environ file in the ../bin dir, default is gfortran"
     echo
     echo "--use_netcdf4  use NetCDF4, the default is NetCDF4. Warning: many of the standard experiments don't work with NetCDF4."
     echo
@@ -87,8 +89,12 @@ if ( $type == EBM ) then
     set cppDefs  = ( "-Duse_netCDF -Duse_netCDF3 -Duse_libMPI -DLAND_BND_TRACERS -DOVERLOAD_C8 -DOVERLOAD_C4 -DOVERLOAD_R4" )
 else if( $type == ACCESS-OM ) then
     set cppDefs  = ( "-Duse_netCDF -Duse_libMPI -DACCESS -DACCESS_OM" )
+else if( $type == ACCESS-OM-BGC ) then
+    set cppDefs  = ( "-Duse_netCDF -Duse_libMPI -DACCESS -DACCESS_OM -DCSIRO_BGC" )
 else if( $type == ACCESS-CM ) then
     set cppDefs  = ( "-Duse_netCDF -Duse_libMPI -DACCESS -DACCESS_CM" )
+else if( $type == ACCESS-ESM ) then
+    set cppDefs  = ( "-Duse_netCDF -Duse_libMPI -DACCESS -DACCESS_CM -DCSIRO_BGC" )
 endif
 
 if ( $unit_testing ) then
@@ -134,7 +140,7 @@ cd $root/exp
 source ./ocean_compile.csh
 if ( $status ) exit $status
 
-if( $type != MOM_solo && $type != ACCESS-OM  && $type != ACCESS-CM) then
+if( $type != MOM_solo && $type != ACCESS-OM  && $type != ACCESS-CM && $type != ACCESS-OM-BGC  && $type != ACCESS-ESM) then
     cd $root/exp
     source ./ice_compile.csh
     if ( $status ) exit $status
@@ -186,12 +192,12 @@ cd $executable:h
 if( $type == MOM_solo ) then
     set srcList = ( mom5/drivers )
     set libs = "$executable:h:h/lib_ocean/lib_ocean.a $executable:h:h/lib_FMS/lib_FMS.a"
-else if( $type == ACCESS-CM ) then
+else if( $type == ACCESS-CM || $type == ACCESS-ESM) then
     set srcList = ( accesscm_coupler )
     set includes = "-I$executable:h:h/lib_FMS -I$executable:h:h/$type/lib_ocean"
     set libs = "$executable:h:h/$type/lib_ocean/lib_ocean.a $executable:h:h/lib_FMS/lib_FMS.a"
     setenv OASIS true
-else if( $type == ACCESS-OM ) then
+else if( $type == ACCESS-OM || $type == ACCESS-OM-BGC) then
     set srcList = ( accessom_coupler )
     set includes = "-I$executable:h:h/lib_FMS -I$executable:h:h/$type/lib_ocean"
     set libs = "$executable:h:h/$type/lib_ocean/lib_ocean.a $executable:h:h/lib_FMS/lib_FMS.a"

exp/ocean_compile.csh

@@ -4,8 +4,11 @@ set srcList = ( mom5/ocean_core mom5/ocean_diag mom5/ocean_wave mom5/ocean_blobs
 
 set lib_name = "lib_ocean"
 
-if( $type == ACCESS-OM || $type == ACCESS-CM ) then
+if( $type == ACCESS-OM || $type == ACCESS-CM || $type == ACCESS-OM-BGC || $type == ACCESS-ESM) then
     set srcList = ( $srcList mom5/ocean_access )
+    if( $type ==  ACCESS-OM-BGC || $type == ACCESS-ESM) then
+        set srcList = ( $srcList mom5/ocean_csiro_bgc )
+    endif
     mkdir -p $executable:h:h/$type/$lib_name
     cd $executable:h:h/$type/$lib_name
 else

src/accessom_coupler/mom_oasis3_interface.F90

@@ -143,11 +143,7 @@ module mom_oasis3_interface_mod
 integer iisc,iiec,jjsc,jjec
 integer iisd,iied,jjsd,jjed
 
-#if defined(ACCESS_WND)
 integer, parameter :: max_fields_in=23
-#else
-integer, parameter :: max_fields_in=22
-#endif
 
 integer, parameter :: max_fields_out=8
 
@@ -360,9 +356,7 @@ subroutine coupler_init(Dom, Time, Time_step_coupled, dt_cpld, Run_len, &
   mom_name_read(20)='wfiform'  !Water flux due to ice forming 
   mom_name_read(21)='licefw'  ! Water flux from land ice
   mom_name_read(22)='liceht'  ! Heat flux from land ice
-#if defined(ACCESS_WND)
   mom_name_read(23)='wnd_io'  !
-#endif
 
   !ocn ==> ice
   mom_name_write(:)=''
@@ -747,10 +741,8 @@ subroutine from_coupler(step,Ocean_sfc,Ice_ocean_boundary, Time)
      Ice_ocean_boundary%licefw(iisc:iiec,jjsc:jjec) =  vwork(iisc:iiec,jjsc:jjec)
   case('liceht')
      Ice_ocean_boundary%liceht(iisc:iiec,jjsc:jjec) =  vwork(iisc:iiec,jjsc:jjec)
-#if defined(ACCESS_WND)
   case('wnd_io')
      Ice_ocean_boundary%wnd(iisc:iiec,jjsc:jjec) =  vwork(iisc:iiec,jjsc:jjec)
-#endif
   case DEFAULT
 ! Probable error. Leave as warning for the moment. RASF
    call mpp_error(WARNING,&

src/mom5/ocean_core/ocean_model.F90

@@ -1638,13 +1638,22 @@ subroutine update_ocean_model(Ice_ocean_boundary, Ocean_state, Ocean_sfc, &
 
        ! compute ocean tendencies from tracer packages
        call mpp_clock_begin(id_otpm_source)
+#if ! defined (CSIRO_BGC)
        call ocean_tpm_source(isd, ied, jsd, jed, Domain, Grid, T_prog(:), T_diag(:), &
             Time, Thickness, Dens, surf_blthick, dtts)
+#else
+       call ocean_tpm_source(isd, ied, jsd, jed, Domain, Grid, T_prog(:), T_diag(:), &
+            Time, Thickness, Dens, surf_blthick, dtts, swflx, sw_frac_zt)
+#endif
        call mpp_clock_end(id_otpm_source)
 
        ! set ocean surface boundary conditions for the tracer packages
        call mpp_clock_begin(id_otpm_bbc)
+#if ! defined (CSIRO_BGC)
        call ocean_tpm_bbc(Domain, Grid, T_prog(:))
+#else
+       call ocean_tpm_bbc(Domain, Grid, T_prog(:), Time)
+#endif
        call mpp_clock_end(id_otpm_bbc)
 
        ! add sponges to T_prog%th_tendency 

src/mom5/ocean_core/ocean_sbc.F90

@@ -534,6 +534,9 @@ module ocean_sbc_mod
 use ocean_util_mod,           only: diagnose_2d, diagnose_2d_u, diagnose_3d_u, diagnose_sum
 use ocean_tracer_util_mod,    only: diagnose_3d_rho
 
+#if defined(CSIRO_BGC)
+use csiro_bgc_mod,            only: csiro_bgc_virtual_fluxes, do_csiro_bgc
+#endif
 implicit none
 
 private
@@ -671,8 +674,11 @@ module ocean_sbc_mod
 #if defined(ACCESS_CM) || defined(ACCESS_OM)
 integer :: id_wfimelt        =-1
 integer :: id_wfiform        =-1
+integer :: id_aice               =-1
+integer :: id_wnd                =-1
 integer :: id_licefw        =-1
 integer :: id_liceht        =-1
+integer :: id_atm_co2            =-1
 #endif
 
 
@@ -815,6 +821,11 @@ module ocean_sbc_mod
 real, allocatable, dimension(:,:,:) :: sslope
 real, allocatable, dimension(:,:) :: aice
 #endif
+#if defined(ACCESS_CM)
+real, allocatable, dimension(:,:) :: co2flux
+real, allocatable, dimension(:,:) :: ocn_co2
+real, allocatable, dimension(:,:) :: atm_co2
+#endif
 
 
 ! ice-ocean-boundary fields are allocated using absolute
@@ -1102,10 +1113,12 @@ subroutine ocean_sbc_init(Grid, Domain, Time, T_prog, T_diag, &
 #if defined(ACCESS_CM) || defined(ACCESS_OM)
   allocate ( Ocean_sfc%gradient (isc_bnd:iec_bnd,jsc_bnd:jec_bnd,2))
   allocate ( sslope(isc:iec, jsc:jec, 2) )
-  allocate ( aice(isc:iec, jsc:jec) )
+  allocate ( aice(isd:ied, jsd:jed) )
 #if defined(ACCESS_CM)
   allocate ( Ocean_sfc%co2    (isc_bnd:iec_bnd,jsc_bnd:jec_bnd), &
              Ocean_sfc%co2flux (isc_bnd:iec_bnd,jsc_bnd:jec_bnd)) 
+  allocate ( co2flux(isd:ied,jsd:jed),ocn_co2(isd:ied,jsd:jed))
+  allocate ( atm_co2(isd:ied,jsd:jed))
 #endif
 #endif
 
@@ -1124,6 +1137,9 @@ subroutine ocean_sbc_init(Grid, Domain, Time, T_prog, T_diag, &
 #if defined(ACCESS_CM)
   Ocean_sfc%co2       = 0.0 
   Ocean_sfc%co2flux   = 0.0 
+  co2flux             = 0.0
+  ocn_co2             = 0.0
+  atm_co2             = 0.0
 #endif
 
   Ocean_sfc%area    = Grid%dat(isc:iec, jsc:jec) * Grid%tmask(isc:iec, jsc:jec, 1) !grid cell area
@@ -1965,12 +1981,26 @@ subroutine ocean_sbc_diag_init(Time, Dens, T_prog)
        '(kg/sec)/1e15', missing_value=missing_value,range=(/-1.e10,1.e10/))   
 
 #if defined(ACCESS_CM) || defined(ACCESS_OM)
- id_total_ocean_wfimelt = register_diag_field('ocean_model','total_ocean_wfimelt',  &
+  id_aice = register_diag_field('ocean_model','aice', Grd%tracer_axes(1:2),&
+       Time%model_time, 'fraction of surface area covered with ice', 'm^2/m^2' ,  &
+       missing_value=missing_value,range=(/-1.e1,1.e1/),                      &
+       standard_name='areal_ice_concentration' )
+  id_wnd = register_diag_field('ocean_model','wnd', Grd%tracer_axes(1:2),&
+       Time%model_time, 'Wind speed', 'm/s' ,  &
+       missing_value=missing_value,range=(/-1.e3,1.e3/),                      &
+       standard_name='wind_speed' )
+  id_total_ocean_wfimelt = register_diag_field('ocean_model','total_ocean_wfimelt',  &
        Time%model_time, 'total icemelt into ocean (>0 enters ocean)',     &
        'kg/sec/1e15', missing_value=missing_value,range=(/-1.e10,1.e10/))
- id_total_ocean_wfiform = register_diag_field('ocean_model','total_ocean_wfiform',  &
+  id_total_ocean_wfiform = register_diag_field('ocean_model','total_ocean_wfiform',  &
        Time%model_time, 'total iceform outof ocean (>0 enters ocean)',     &
        'kg/sec/1e15', missing_value=missing_value,range=(/-1.e10,1.e10/))
+#if defined(ACCESS_CM)
+  id_atm_co2 = register_diag_field('ocean_model','atm_co2', Grd%tracer_axes(1:2),&
+       Time%model_time, 'Atmospheric CO2 content', 'ppm' ,  &
+       missing_value=missing_value,range=(/-1.e1,1.e4/),                      &
+       standard_name='atmospheric_co2' )
+#endif
   id_total_ocean_licefw = register_diag_field('ocean_model','total_ocean_licefw',  &
         Time%model_time, 'total land icemelt into ocean (>0 enters ocean)',     &
         'kg/sec/1e15', missing_value=missing_value,range=(/-1.e10,1.e10/))
@@ -2773,6 +2803,10 @@ subroutine initialize_ocean_sfc(Time, Thickness, T_prog, T_diag, Velocity, Ocean
       Ocean_sfc%v_surf(isc_bnd:iec_bnd,jsc_bnd:jec_bnd)  = Velocity%u(isc:iec,jsc:jec,1,2,taup1)
       Ocean_sfc%sea_lev(isc_bnd:iec_bnd,jsc_bnd:jec_bnd) = Thickness%sea_lev(isc:iec,jsc:jec)
       Ocean_sfc%frazil(isc_bnd:iec_bnd,jsc_bnd:jec_bnd)  = 0.0
+#if defined(ACCESS_CM)
+      Ocean_sfc%co2flux(isc_bnd:iec_bnd,jsc_bnd:jec_bnd)     = co2flux(isc:iec,jsc:jec)  !These should really come from a restart RASF
+      Ocean_sfc%co2(isc_bnd:iec_bnd,jsc_bnd:jec_bnd)         = ocn_co2(isc:iec,jsc:jec)
+#endif
   end where
 
   ! when enabled, use FAFMIP redistributed heat tracer for sst
@@ -2807,6 +2841,11 @@ subroutine initialize_ocean_sfc(Time, Thickness, T_prog, T_diag, Velocity, Ocean
   id_field = register_restart_field(Sfc_restart, filename, 'v_surf', Ocean_sfc%v_surf,Ocean_sfc%Domain)
   id_field = register_restart_field(Sfc_restart, filename, 'sea_lev',Ocean_sfc%sea_lev,Ocean_sfc%Domain)
   id_field = register_restart_field(Sfc_restart, filename, 'frazil', Ocean_sfc%frazil,Ocean_sfc%Domain)
+#if defined(ACCESS_CM)
+!RASF Make these optional so we don't break existing runs.
+  id_field = register_restart_field(Sfc_restart, filename, 'co2flux',Ocean_sfc%co2flux,Ocean_sfc%Domain, mandatory=.false.)
+  id_field = register_restart_field(Sfc_restart, filename, 'ocn_co2', Ocean_sfc%co2,Ocean_sfc%Domain, mandatory=.false.)
+#endif
   if (file_exist('INPUT/ocean_sbc.res.nc')) then
      call restore_state(Sfc_restart)
   endif
@@ -2891,6 +2930,10 @@ subroutine sum_ocean_sfc(Time, Thickness, T_prog, T_diag, Dens, Velocity, Ocean_
             Ocean_sfc%sea_lev(i,j) = Ocean_sfc%sea_lev(i,j) + Thickness%sea_lev(ii,jj)
 #if defined(ACCESS_CM) || defined(ACCESS_OM)
             Ocean_sfc%gradient(i,j,:) = Ocean_sfc%gradient(i,j,:) + sslope(ii,jj,:)  
+#endif
+#if defined(ACCESS_CM)
+            Ocean_sfc%co2flux(i,j) = Ocean_sfc%co2flux(i,j) + co2flux(ii,jj)
+            Ocean_sfc%co2(i,j)     = Ocean_sfc%co2(i,j) + ocn_co2(ii,jj)
 #endif
          enddo
       enddo
@@ -2990,6 +3033,10 @@ subroutine zero_ocean_sfc(Ocean_sfc)
         Ocean_sfc%frazil(i,j) = 0.0
 #if defined(ACCESS_CM) || defined(ACCESS_OM)
         Ocean_sfc%gradient(i,j,:)= 0.0
+#endif
+#if defined(ACCESS_CM)
+        Ocean_sfc%co2flux(i,j) = 0.0
+        Ocean_sfc%co2(i,j)     = 0.0
 #endif
      enddo
   enddo
@@ -3067,6 +3114,10 @@ subroutine avg_ocean_sfc(Time, Thickness, T_prog, T_diag, Velocity, Ocean_sfc)
            Ocean_sfc%sea_lev(i,j) = Ocean_sfc%sea_lev(i,j)*divid 
 #if defined(ACCESS_CM) || defined(ACCESS_OM)
            Ocean_sfc%gradient(i,j,:) = Ocean_sfc%gradient(i,j,:)*divid
+#endif
+#if defined(ACCESS_CM)
+           Ocean_sfc%co2flux(i,j) = Ocean_sfc%co2flux(i,j)*divid
+           Ocean_sfc%co2(i,j)     = Ocean_sfc%co2(i,j)*divid
 #endif
         endif
      enddo
@@ -3092,6 +3143,10 @@ subroutine avg_ocean_sfc(Time, Thickness, T_prog, T_diag, Velocity, Ocean_sfc)
             if(Grd%tmask(ii,jj,1) == 1.0) then
                 Ocean_sfc%s_surf(i,j) = T_prog(index_salt)%field(ii,jj,1,taup1)
                 Ocean_sfc%sea_lev(i,j)= Thickness%sea_lev(ii,jj)
+#if defined(ACCESS_CM)
+                Ocean_sfc%co2flux(i,j) = co2flux(ii,jj)
+                Ocean_sfc%co2(i,j)     = ocn_co2(ii,jj)
+#endif
             endif
          enddo
       enddo
@@ -4270,13 +4325,43 @@ subroutine get_ocean_sbc(Time, Ice_ocean_boundary, Thickness, Dens, Ext_mode, T_
      call mpp_update_domains(pme(:,:), Dom%domain2d)
   endif 
 
+#if defined(ACCESS_CM)
+  do j = jsc_bnd,jec_bnd
+     do i = isc_bnd,iec_bnd
+        ii = i + i_shift
+        jj = j + j_shift
+        atm_co2(ii,jj) = Ice_ocean_boundary%co2(i,j)*Grd%tmask(ii,jj,1)
+     enddo
+  enddo
+#endif
+
+#if defined(ACCESS_CM) || defined(ACCESS_OM)
+  do j = jsc_bnd,jec_bnd
+     do i = isc_bnd,iec_bnd
+        ii = i + i_shift
+        jj = j + j_shift
+        aice(ii,jj) = Ice_ocean_boundary%aice(i,j)*Grd%tmask(ii,jj,1)
+     enddo
+  enddo
+#endif
 
 
   !--------compute surface tracer fluxes from tracer packages------------------- 
   !
+#if defined(ACCESS_CM) && defined(CSIRO_BGC)
+  call ocean_tpm_sbc(Dom, Grd, T_prog(:), Time, Ice_ocean_boundary%fluxes, runoff, &
+                     isc_bnd, iec_bnd, jsc_bnd, jec_bnd,aice, Velocity%u10, &
+     use_waterflux, salt_restore_as_salt_flux, atm_co2, co2flux, ocn_co2)
+#elif defined(ACCESS_OM) && defined(CSIRO_BGC)
+! Do not pass co2flux, ocn_co2 or atm_co2
+  call ocean_tpm_sbc(Dom, Grd, T_prog(:), Time, Ice_ocean_boundary%fluxes, runoff, &
+                     isc_bnd, iec_bnd, jsc_bnd, jec_bnd,aice=aice, wnd=Velocity%u10, &
+                     use_waterflux=use_waterflux, salt_restore_as_salt_flux=salt_restore_as_salt_flux)
+#else 
   call ocean_tpm_sbc(Dom, Grd, T_prog(:), Time, Ice_ocean_boundary%fluxes, runoff, &
                      isc_bnd, iec_bnd, jsc_bnd, jec_bnd)
-
+! Leave this case blank at the moment. We want to use the bgc with SIS etc.
+#endif
 
   !--------send diagnostics------------------- 
   !
@@ -4285,16 +4370,6 @@ subroutine get_ocean_sbc(Time, Ice_ocean_boundary, Thickness, Dens, Ext_mode, T_
                       melt, liquid_precip, frozen_precip, evaporation, sensible, longwave, &
                       latent, swflx, swflx_vis)
 
-#if defined(ACCESS_CM) || defined(ACCESS_OM)
-  do j = jsc_bnd,jec_bnd
-     do i = isc_bnd,iec_bnd
-        ii = i + i_shift
-        jj = j + j_shift
-        aice(ii,jj) = Ice_ocean_boundary%aice(i,j)*Grd%tmask(ii,jj,1)
-     enddo
-  enddo
-#endif
-
 end subroutine get_ocean_sbc
 ! </SUBROUTINE> NAME="get_ocean_sbc"
 
@@ -4536,6 +4611,14 @@ subroutine flux_adjust(Time, T_diag, Dens, Ext_mode, T_prog, Velocity, river, me
              enddo
           enddo
 
+#if defined(CSIRO_BGC)
+          ! if salt fluxes are used to restore salinity, then virtual fluxes are
+          ! needed for csiro BGC tracers. mac, dec12.
+          if (do_csiro_bgc) then
+            call csiro_bgc_virtual_fluxes(isc, iec, jsc, jec, isd, ied, jsd, jed, flx_restore, T_prog)
+          endif
+#endif
+
       endif  ! endif for if (use_waterflux .and. .not. salt_restore_as_salt_flux) then
 
   endif  ! endif for if(id_restore(index_salt) > 0 .and. salt_damp_factor > 0.0)
@@ -5875,7 +5958,12 @@ subroutine ocean_sbc_diag(Time, Velocity, Thickness, Dens, T_prog, Ice_ocean_bou
       total_stuff  = mpp_global_sum(Dom%domain2d,wrk1_2d(:,:), NON_BITWISE_EXACT_SUM)
       used = send_data (id_total_ocean_wfiform, total_stuff*1e-15, Time%model_time)
   endif
-
+  if (id_aice > 0) used = send_data(id_aice, aice(:,:),    &
+                 Time%model_time, rmask=Grd%tmask(:,:,1),  &
+                 is_in=isc, js_in=jsc, ie_in=iec, je_in=jec)
+  if (id_wnd > 0) used = send_data(id_wnd, Velocity%u10(:,:),    &
+                 Time%model_time, rmask=Grd%tmask(:,:,1),  &
+                 is_in=isc, js_in=jsc, ie_in=iec, je_in=jec)
    ! waterflux associated with land ice melt into ocean (kg/(m2*sec))
    if (id_licefw > 0) then
        do j=jsc_bnd,jec_bnd
@@ -5889,6 +5977,11 @@ subroutine ocean_sbc_diag(Time, Velocity, Thickness, Dens, T_prog, Ice_ocean_bou
               Time%model_time, rmask=Grd%tmask(:,:,1),  &
               is_in=isc, js_in=jsc, ie_in=iec, je_in=jec)
    endif
+#if defined(ACCESS_CM)
+  if (id_atm_co2 > 0) used = send_data(id_atm_co2, atm_co2(:,:),    &
+                 Time%model_time, rmask=Grd%tmask(:,:,1),  &
+                 is_in=isc, js_in=jsc, ie_in=iec, je_in=jec)
+#endif
 #endif