<p><b>duda</b> 2010-04-28 17:25:27 -0600 (Wed, 28 Apr 2010)</p><p>BRANCH COMMIT<br>
<br>
Add implementation of mpas_restart_to_cam() routine plus code<br>
in the test driver to exercise the restart capabilities of the<br>
interface.<br>
<br>
<br>
M module_mpas_cam_interface.F<br>
M test_driver.F<br>
</p><hr noshade><pre><font color="gray">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-04-28 23:11:14 UTC (rev 221)
+++ branches/mpas_cam_coupling/src/driver_cam_interface/module_mpas_cam_interface.F 2010-04-28 23:25:27 UTC (rev 222)
@@ -373,8 +373,96 @@
block % time_levs(1) % state % xtime % scalar = 0.0
end subroutine cam_inidat_to_mpas
+
+
+ subroutine mpas_restart_to_cam(Numcols, MaxEdges, Plev, Pcnst, Psd, Phis, Theta, U, W, Tracer)
+ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ ! SUBROUTINE MPAS_RESTART_TO_CAM
+ !
+ ! Input: Numcols - number of columns/cells
+ ! MaxEdges - maximum number of edges per cell
+ ! Plev - number of vertical levels
+ ! Pcnst - number of tracers
+ ! Psd - dry surface pressure (Pa)
+ ! Phis - surface geopotential (m^2/s^2)
+ ! Theta - potential temperature (K)
+ ! U - normal velocity at edges (m/s)
+ ! W - vertical velocity (Pa/s)
+ ! Tracer - tracer mixing ratios (kg/kg)
+ !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+ implicit none
+
+ integer, intent(in) :: Numcols
+ integer, intent(in) :: MaxEdges
+ integer, intent(in) :: Plev
+ integer, intent(in) :: Pcnst
+ real (kind=RKIND), dimension(Numcols), intent(out) :: Psd
+ real (kind=RKIND), dimension(Numcols), intent(out) :: Phis
+ real (kind=RKIND), dimension(Numcols,Plev), intent(out) :: Theta
+ real (kind=RKIND), dimension(Numcols,MaxEdges,Plev), intent(out) :: U
+ real (kind=RKIND), dimension(Numcols,Plev+1), intent(out) :: W
+ real (kind=RKIND), dimension(Numcols,Plev,Pcnst), intent(out) :: Tracer
+
+ integer :: iCell, iEdge, k, iScalar, i
+ type (block_type), pointer :: block
+ integer, dimension(:), pointer :: nEdgesOnCell
+ integer, dimension(:,:), pointer :: edgesOnCell
+ write(0,*) 'Called MPAS_RESTART_TO_CAM'
+
+ block => domain % blocklist
+
+ nEdgesOnCell => block % mesh % nEdgesOnCell % array
+ edgesOnCell => block % mesh % edgesOnCell % array
+
+ !
+ ! Perform basic sanity check on expected and available field dimensions
+ !
+ if (Numcols /= block % mesh % nCellsSolve) then
+ write(0,*) 'Error: mismatch between Numcols and nCellsSolve: ', Numcols, block % mesh % nCellsSolve
+ return
+ end if
+ if (Plev /= block % mesh % nVertLevels) then
+ write(0,*) 'Error: mismatch between Plev and nVertLevels: ', Plev, block % mesh % nVertLevels
+ return
+ end if
+ if (Pcnst /= num_scalars) then
+ write(0,*) 'Error: mismatch between Pcnst and num_scalars: ', Pcnst, num_scalars
+ return
+ end if
+
+
+ !
+ ! Fill in output arrays from block % time_levs(1) % state arrays
+ ! NOTE: We must copy from time_levs(1) rather than time_levs(2), since
+ ! mpas_dyn_run swaps the new time level into time_levs(1) after
+ ! taking a timestep.
+ !
+ do iCell=1,block % mesh % nCellsSolve
+ Psd(iCell) = block % time_levs(1) % state % surface_pressure % array(iCell)
+ Phis(iCell) = block % time_levs(1) % state % geopotential % array(1,iCell)
+ do k=1,block % mesh % nVertLevels
+ Theta(iCell,k) = block % time_levs(1) % state % theta % array(k,iCell)
+ W(iCell,k) = block % time_levs(1) % state % ww % array(k,iCell)
+ !
+ ! NOTE: Eventually, we need to ensure that the moisture variables we return to CAM
+ ! are what CAM expects, rather than simply what we have in our scalars array
+ !
+ do iScalar=1,num_scalars
+ Tracer(iCell,k,iScalar) = block % time_levs(1) % state % scalars % array(iScalar,k,iCell)
+ end do
+ do i=1,nEdgesOnCell(iCell)
+ iEdge = edgesOnCell(i,iCell)
+ U(iCell,i,k) = block % time_levs(1) % state % u % array(k,iEdge)
+ end do
+ end do
+ W(iCell,k) = block % time_levs(1) % state % ww % array(k,iCell)
+ end do
+
+ end subroutine mpas_restart_to_cam
+
+
subroutine cam_restart_to_mpas(Numcols, MaxEdges, Plev, Pcnst, Pref, Ptop, &
Ae, Be, Ac, Bc, Psd, Phis, Theta, U, W, Tracer)
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
Modified: branches/mpas_cam_coupling/src/driver_cam_interface/test_driver.F
===================================================================
--- branches/mpas_cam_coupling/src/driver_cam_interface/test_driver.F 2010-04-28 23:11:14 UTC (rev 221)
+++ branches/mpas_cam_coupling/src/driver_cam_interface/test_driver.F 2010-04-28 23:25:27 UTC (rev 222)
@@ -13,9 +13,12 @@
real (kind=RKIND), dimension(26) :: Ac, Bc
real (kind=RKIND), dimension(27) :: Ae, Be
real (kind=RKIND), dimension(:), pointer :: Psd, Phis
+ real (kind=RKIND), dimension(:), pointer :: Psd_rst, Phis_rst
real (kind=RKIND), dimension(:,:), pointer :: T, Ux, Uy, dPdt
real (kind=RKIND), dimension(:,:), pointer :: T_tend, Ux_tend, Uy_tend
+ real (kind=RKIND), dimension(:,:), pointer :: Theta_rst, W_rst
real (kind=RKIND), dimension(:,:,:), pointer :: U, Tracer
+ real (kind=RKIND), dimension(:,:,:), pointer :: U_rst, Tracer_rst
call MPI_Init(mpi_ierr)
call MPI_Comm_rank(MPI_COMM_WORLD, mpi_procID, mpi_ierr)
@@ -71,6 +74,13 @@
allocate(Ux_tend(nCellsSolve,nVertLevels))
allocate(Uy_tend(nCellsSolve,nVertLevels))
+ allocate(Psd_rst(nCellsSolve))
+ allocate(Phis_rst(nCellsSolve))
+ allocate(Theta_rst(nCellsSolve,nVertLevels))
+ allocate(U_rst(nCellsSolve,maxEdges,nVertLevels))
+ allocate(W_rst(nCellsSolve,nVertLevels+1))
+ allocate(Tracer_rst(nCellsSolve,nVertLevels,nScalars))
+
read(22) Psd
read(22) Phis
read(22) T
@@ -86,19 +96,41 @@
! Begin calling CAM-MAS interface routines
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
+
+ !
+ ! Make "cold-start" run
+ !
call mpas_init1(MPI_COMM_WORLD)
call cam_inidat_to_mpas(nCellsSolve, maxEdges, nVertLevels, nScalars, Pref, Ptop, &
Ae, Be, Ac, Bc, Psd, Phis, T, U, Tracer)
+
call mpas_init2()
- do istep=1,24
+ do istep=1,12
call mpas_dyn_run()
end do
+ call mpas_restart_to_cam(nCellsSolve, maxEdges, nVertLevels, nScalars, &
+ Psd_rst, Phis_rst, Theta_rst, U_rst, W_rst, Tracer_rst)
+ write(filename, '(a,i4.4)') 'restart.', mpi_procID
+ open(22, file=trim(filename), form='unformatted', status='unknown')
+ write(22) Psd_rst
+ write(22) Phis_rst
+ write(22) Theta_rst
+ write(22) U_rst
+ write(22) W_rst
+ write(22) Tracer_rst
+ close(22)
+
+ do istep=1,12
+ call mpas_dyn_run()
+ end do
+
! *** Test routine to write model state -- not part of MPAS-CAM interface
call mpas_write_output()
+ if (mpi_procID == 0) call system('mv output.nc output1.nc')
call mpas_to_cam(nCellsSolve, nVertLevels, nScalars, &
Psd, Phis, T, Ux, Uy, dPdt, Tracer)
@@ -108,6 +140,38 @@
call mpas_final()
+
+ !
+ ! Make "restart" run
+ !
+ call mpas_init1(MPI_COMM_WORLD)
+
+ write(filename, '(a,i4.4)') 'restart.', mpi_procID
+ open(22, file=trim(filename), form='unformatted', status='old')
+ read(22) Psd_rst
+ read(22) Phis_rst
+ read(22) Theta_rst
+ read(22) U_rst
+ read(22) W_rst
+ read(22) Tracer_rst
+ close(22)
+
+ call cam_restart_to_mpas(nCellsSolve, maxEdges, nVertLevels, nScalars, Pref, Ptop, &
+ Ae, Be, Ac, Bc, Psd_rst, Phis_rst, Theta_rst, U_rst, W_rst, Tracer_rst)
+
+ call mpas_init2()
+
+ do istep=1,12
+ call mpas_dyn_run()
+ end do
+
+ ! *** Test routine to write model state -- not part of MPAS-CAM interface
+ call mpas_write_output()
+ if (mpi_procID == 0) call system('mv output.nc output2.nc')
+
+ call mpas_final()
+
+
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
! Finish calling CAM-MAS interface routines
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
@@ -125,6 +189,13 @@
deallocate(Ux_tend)
deallocate(Uy_tend)
+ deallocate(Psd_rst)
+ deallocate(Phis_rst)
+ deallocate(Theta_rst)
+ deallocate(U_rst)
+ deallocate(W_rst)
+ deallocate(Tracer_rst)
+
call MPI_Finalize(mpi_ierr)
stop
</font>
</pre>