<p><b>duda</b> 2010-05-04 14:26:01 -0600 (Tue, 04 May 2010)</p><p>BRANCH COMMIT<br>
<br>
Update interface for MPAS_INIT1 to accept a second argument, which is<br>
the physics timestep, and compute and store the scalar tendendies, rather<br>
than storing the new scalar fields returned by CAM physics.<br>
<br>
<br>
M src/core_hyd_atmos/Registry<br>
M src/driver_cam_interface/module_mpas_cam_interface.F<br>
M src/driver_cam_interface/test_driver.F<br>
</p><hr noshade><pre><font color="gray">Modified: branches/mpas_cam_coupling/src/core_hyd_atmos/Registry
===================================================================
--- branches/mpas_cam_coupling/src/core_hyd_atmos/Registry        2010-05-04 20:25:38 UTC (rev 242)
+++ branches/mpas_cam_coupling/src/core_hyd_atmos/Registry        2010-05-04 20:26:01 UTC (rev 243)
@@ -155,6 +155,9 @@
# Physics tendencies
var real theta_phys_tend ( nVertLevels nCells Time ) - theta_phys_tend - -
var real u_phys_tend ( nVertLevels nEdges Time ) - u_phys_tend - -
+var real qv_phys_tend ( nVertLevels nCells Time ) - qv_phys_tend scalars_phys_tend moist_phys_tend
+var real qc_phys_tend ( nVertLevels nCells Time ) - qc_phys_tend scalars_phys_tend moist_phys_tend
+var real qr_phys_tend ( nVertLevels nCells Time ) - qr_phys_tend scalars_phys_tend moist_phys_tend
# Space needed for advection
var real deriv_two ( FIFTEEN TWO nEdges ) o deriv_two - -
Modified: branches/mpas_cam_coupling/src/driver_cam_interface/module_mpas_cam_interface.F
===================================================================
--- branches/mpas_cam_coupling/src/driver_cam_interface/module_mpas_cam_interface.F        2010-05-04 20:25:38 UTC (rev 242)
+++ branches/mpas_cam_coupling/src/driver_cam_interface/module_mpas_cam_interface.F        2010-05-04 20:26:01 UTC (rev 243)
@@ -10,14 +10,14 @@
type (dm_info), pointer :: dminfo
type (domain_type), pointer :: domain
- real (kind=RKIND) :: dt
+ real (kind=RKIND) :: dt, dt_physics
integer :: itimestep
contains
- subroutine mpas_init1(mpi_comm)
+ subroutine mpas_init1(mpi_comm, phys_dt)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! SUBROUTINE MPAS_INIT1
!
@@ -33,6 +33,7 @@
implicit none
integer, intent(in) :: mpi_comm
+ real (kind=RKIND), intent(in) :: phys_dt
integer :: iCell, iEdge, j
real (kind=8) :: m, angle
@@ -49,6 +50,16 @@
call read_namelist(dminfo)
+ !
+ ! Set physics timestep and verify that it is evenly divided by dynamics timestep
+ !
+ dt_physics = phys_dt
+ if (mod(dt_physics, config_dt) /= 0.) then
+ write(0,*) 'Error: physics timestep does not evenly divide dynamics timestep'
+ return
+ end if
+
+
call allocate_domain(domain, dminfo)
! NOTE: Really, we only want to read in the grid information, since any state
@@ -713,7 +724,7 @@
integer, dimension(:,:), pointer :: edgesOnCell
real (kind=RKIND), dimension(:,:), pointer :: east, north, edge_normal, pressure
real (kind=RKIND), dimension(:,:), pointer :: Ux_tend_halo, Uy_tend_halo, theta_tend, u_tend
- real (kind=RKIND), dimension(:,:,:), pointer :: scalars
+ real (kind=RKIND), dimension(:,:,:), pointer :: scalars, scalars_phys_tend
type (block_type), pointer :: block
write(0,*) 'Called CAM_TO_MPAS'
@@ -726,6 +737,7 @@
edge_normal => block % mesh % edge_normal % array
pressure => block % time_levs(2) % state % pressure % array
scalars => block % time_levs(2) % state % scalars % array
+ scalars_phys_tend => block % time_levs(2) % state % scalars_phys_tend % array
theta_tend => block % time_levs(2) % state % theta_phys_tend % array
u_tend => block % time_levs(2) % state % u_phys_tend % array
@@ -793,7 +805,7 @@
! properly assign variables from Tracer array provided from CAM physics
!
do iScalar=1,num_scalars
- scalars(iScalar,k,iCell) = Tracer(iCell,k,iScalar)
+ scalars_phys_tend(iScalar,k,iCell) = (Tracer(iCell,k,iScalar) - scalars(iScalar,k,iCell)) / dt_physics
end do
end do
end do
@@ -832,7 +844,7 @@
! Move time level 2 fields back into time level 1 for next time step
call shift_time_levels(domain)
-
+
end subroutine mpas_dyn_run
Modified: branches/mpas_cam_coupling/src/driver_cam_interface/test_driver.F
===================================================================
--- branches/mpas_cam_coupling/src/driver_cam_interface/test_driver.F        2010-05-04 20:25:38 UTC (rev 242)
+++ branches/mpas_cam_coupling/src/driver_cam_interface/test_driver.F        2010-05-04 20:26:01 UTC (rev 243)
@@ -9,6 +9,7 @@
character (len=32) :: filename
integer :: nCellsSolve, maxEdges, nVertLevels, nScalars, istep
+ real (kind=RKIND) :: phys_dt
real (kind=RKIND) :: Pref, Ptop
real (kind=RKIND), dimension(26) :: Ac, Bc
real (kind=RKIND), dimension(27) :: Ae, Be
@@ -46,6 +47,8 @@
0.620120200, 0.723535500, 0.817676800, 0.896215300, &
0.953476103, 0.985112200, 1.000000000 /
+ phys_dt = 10800.
+
do k=1,26
Ac(k) = 0.5*(Ae(k) + Ae(k+1))
Bc(k) = 0.5*(Be(k) + Be(k+1))
@@ -100,7 +103,7 @@
!
! Make "cold-start" run
!
- call mpas_init1(MPI_COMM_WORLD)
+ call mpas_init1(MPI_COMM_WORLD, phys_dt)
call cam_inidat_to_mpas(nCellsSolve, maxEdges, nVertLevels, nScalars, Pref, Ptop, &
Ae, Be, Ac, Bc, Psd, Phis, T, U, Tracer)
@@ -144,7 +147,7 @@
!
! Make "restart" run
!
- call mpas_init1(MPI_COMM_WORLD)
+ call mpas_init1(MPI_COMM_WORLD, phys_dt)
write(filename, '(a,i4.4)') 'restart.', mpi_procID
open(22, file=trim(filename), form='unformatted', status='old')
</font>
</pre>