<p><b>dwj07@fsu.edu</b> 2013-01-10 14:06:14 -0700 (Thu, 10 Jan 2013)</p><p><br>
        -- BRANCH COMMIT --<br>
<br>
        Adding coupled time averager.<br>
</p><hr noshade><pre><font color="gray">Added: branches/ocean_projects/generic_forcing/src/core_ocean/mpas_ocn_time_average_coupled.F
===================================================================
--- branches/ocean_projects/generic_forcing/src/core_ocean/mpas_ocn_time_average_coupled.F         (rev 0)
+++ branches/ocean_projects/generic_forcing/src/core_ocean/mpas_ocn_time_average_coupled.F        2013-01-10 21:06:14 UTC (rev 2381)
@@ -0,0 +1,119 @@
+module ocn_time_average_coupled
+
+ use mpas_grid_types
+ use ocn_constants
+
+ implicit none
+ save
+ public
+
+ contains
+
+ subroutine ocn_time_average_coupled_init(state)!{{{
+ type (state_type), intent(inout) :: state
+
+ real (kind=RKIND), pointer :: nAccumulateCoupled
+
+ real (kind=RKIND), dimension(:), pointer :: temperatureState, salinityState, meridionalVelocityState, zonalVelocitystate, meridionalSSHGradientState, zonalSSHGradientState
+
+#ifdef MPAS_CESM
+ nAccumulateCoupled => state % nAccumulateCoupled % scalar
+
+ temperatureState => state % temperatureState % array
+ salinityState => state % salinityState % array
+ meridionalVelocityState => state % meridionalVelocityState % array
+ zonalVelocityState => state % zonalVelocityState % array
+ meridionalSSHGradientState => state % meridionalSSHGradientState % array
+ zonalSSHGradientState => state % zonalSSHGradientState % array
+
+ nAccumulateCoupled = 0
+
+ temperatureState = 0.0_RKIND
+ salinityState = 0.0_RKIND
+ meridionalVelocityState = 0.0_RKIND
+ zonalVelocityState = 0.0_RKIND
+ meridionalSSHGradientState = 0.0_RKIND
+ zonalSSHGradientState = 0.0_RKIND
+#endif
+
+ end subroutine ocn_time_average_coupled_init!}}}
+
+ subroutine ocn_time_average_coupled_accumulate(state, old_state)!{{{
+ type (state_type), intent(inout) :: state
+ type (state_type), intent(in) :: old_state
+
+ integer :: index_temperature, index_salinity
+
+ real (kind=RKIND), pointer :: nAccumulateCoupled, old_nAccumulateCoupled
+
+ real (kind=RKIND), dimension(:), pointer :: temperatureState, salinityState, meridionalVelocityState, zonalVelocitystate, meridionalSSHGradientState, zonalSSHGradientState
+ real (kind=RKIND), dimension(:), pointer :: old_temperatureState, old_salinityState, old_meridionalVelocityState, old_zonalVelocitystate, old_meridionalSSHGradientState, old_zonalSSHGradientState
+ real (kind=RKIND), dimension(:), pointer :: ssh
+ real (kind=RKIND), dimension(:,:), pointer :: uReconstructZonal, uReconstructMeridional
+ real (kind=RKIND), dimension(:,:,:), pointer :: tracers
+
+#ifdef MPAS_CESM
+ index_temperature = state % index_temperature
+ index_salinity = state % index_salinity
+
+ nAccumulateCoupled => state % nAccumulateCoupled % scalar
+ old_nAccumulateCoupled => old_state % nAccumulateCoupled % scalar
+
+ nAccumulateCoupled = old_nAccumulateCoupled + 1
+
+ ssh => state % ssh % array
+ uReconstructZonal => state % uReconstructZonal % array
+ uReconstructMeridional => state % uReconstructMeridional % array
+ tracers => state % tracers % array
+
+ temperatureState => state % temperatureState % array
+ salinityState => state % salinityState % array
+ meridionalVelocityState => state % meridionalVelocityState % array
+ zonalVelocityState => state % zonalVelocityState % array
+ meridionalSSHGradientState => state % meridionalSSHGradientState % array
+ zonalSSHGradientState => state % zonalSSHGradientState % array
+
+ old_temperatureState => old_state % temperatureState % array
+ old_salinityState => old_state % salinityState % array
+ old_meridionalVelocityState => old_state % meridionalVelocityState % array
+ old_zonalVelocityState => old_state % zonalVelocityState % array
+ old_meridionalSSHGradientState => old_state % meridionalSSHGradientState % array
+ old_zonalSSHGradientState => old_state % zonalSSHGradientState % array
+
+ temperatureState = old_temperatureState + (tracers(index_temperature, 1, :) + T0_Kelvin)
+ salinityState = old_salinityState + tracers(index_salinity, 1, :)
+ meridionalVelocityState = old_meridionalVelocityState + uReconstructMeridional(1,:)
+ zonalVelocityState = old_zonalVelocityState + uReconstructZonal(1,:)
+#endif
+
+ end subroutine ocn_time_average_coupled_accumulate!}}}
+
+ subroutine ocn_time_average_coupled_normalize(state)!{{{
+ type (state_type), intent(inout) :: state
+
+ real (kind=RKIND), pointer :: nAccumulateCoupled
+
+ real (kind=RKIND), dimension(:), pointer :: temperatureState, salinityState, meridionalVelocityState, zonalVelocitystate, meridionalSSHGradientState, zonalSSHGradientState
+
+#ifdef MPAS_CESM
+ nAccumulateCoupled => state % nAccumulateCoupled % scalar
+
+ temperatureState => state % temperatureState % array
+ salinityState => state % salinityState % array
+ meridionalVelocityState => state % meridionalVelocityState % array
+ zonalVelocityState => state % zonalVelocityState % array
+ meridionalSSHGradientState => state % meridionalSSHGradientState % array
+ zonalSSHGradientState => state % zonalSSHGradientState % array
+
+ if(nAccumulateCoupled > 0) then
+ temperatureState = temperatureState / nAccumulateCoupled
+ salinityState = salinityState / nAccumulateCoupled
+ meridionalVelocityState = meridionalVelocityState / nAccumulateCoupled
+ zonalVelocityState = zonalVelocityState / nAccumulateCoupled
+ meridionalSSHGradientState = meridionalSSHGradientState / nAccumulateCoupled
+ zonalSSHGradientState = zonalSSHGradientState / nAccumulateCoupled
+ end if
+#endif
+ end subroutine ocn_time_average_coupled_normalize!}}}
+
+end module ocn_time_average_coupled
</font>
</pre>