<p><b>mhoffman@lanl.gov</b> 2012-04-18 16:42:58 -0600 (Wed, 18 Apr 2012)</p><p>BRANCH COMMIT<br>
Adjusted Forward Euler time stepper to make velocity a diagnostic field.<br>
Added a diagnostic solve for velocity in init to define an initial velocity condition.<br>
</p><hr noshade><pre><font color="gray">Modified: branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_mpas_core.F
===================================================================
--- branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_mpas_core.F        2012-04-18 21:38:59 UTC (rev 1795)
+++ branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_mpas_core.F        2012-04-18 22:42:58 UTC (rev 1796)
@@ -21,7 +21,7 @@
use mpas_configure
use mpas_grid_types
use mpas_constants
- use land_ice_vel, only: land_ice_vel_init
+ use land_ice_vel
implicit none
@@ -52,8 +52,16 @@
block => domain % blocklist
do while (associated(block))
call mpas_init_block(block, block % mesh, dt)
+ ! Assign initial time stamp
block % state % time_levs(1) % state % xtime % scalar = startTimeStamp
- ! Make sure all time levels have a copy of the initial state prior the first velocity solve.
+
+ ! Calculate an initial (diagnostic) velocity field consistent with the initial thickness field
+ ! \todo: skip this step if a velocity field was supplied in the I.C. input file
+ call land_ice_vel_solve(block % mesh, block % state % time_levs(1) % state, &
+ block % state % time_levs(1) % state % normalVelocity % array, err)
+ ! (halo update of velocity field)
+
+ ! Make sure all time levels have a copy of the initial state
do i=2,nTimeLevs
call mpas_copy_state(block % state % time_levs(i) % state, block % state % time_levs(1) % state)
end do
@@ -63,6 +71,7 @@
current_outfile_frames = 0
+
if(err == 1) then
print *, "An error has occurred in mpas_core_init. Aborting..."
call mpas_dmpar_global_abort("An error has occurred in mpas_core_init. Aborting...")
Modified: branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_time_integration_forwardeuler.F
===================================================================
--- branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_time_integration_forwardeuler.F        2012-04-18 21:38:59 UTC (rev 1795)
+++ branches/land_ice_projects/implement_core/src/core_land_ice/mpas_land_ice_time_integration_forwardeuler.F        2012-04-18 22:42:58 UTC (rev 1796)
@@ -131,20 +131,24 @@
end where
print *,'Cells with nonzero thickness:', sum(mask)
- ! Solve the velocity (using old state)
- call land_ice_vel_solve(mesh, stateOld, normalVelocityNew, err)
- ! update halos on velocity
- ! call ()
+
! === Calculate diagnostic variables for new state =====================
! Remap layer thicknesses (could be moved to diagnostic solve)
do k = 1, nVertLevels
layerThicknessNew(k, :) = thicknessNew * layerThicknessFractions(k)
end do
- ! call land_ice_diagnostic_solve()
+ ! call land_ice_diagnostic_solve() ! perhaps layer thickness remapping and velocity solve should move in here.
- ! Reconstruct velocities for this time step (could be moved to diagnostic solve)
+ ! Compute the diagnostic velocity for this time step using the updated (new) state
+ ! Pass in the old velocity to be used as an initial guess for velocity solvers that need one.
+ normalVelocityNew = normalVelocityOld
+ call land_ice_vel_solve(mesh, stateNew, normalVelocityNew, err)
+ ! update halos on velocity
+ ! call ()
+
+ ! Calculate reconstructed velocities for this time step (could be moved to diagnostic solve or vel_solve)
call mpas_reconstruct(mesh, normalVelocityNew, &
uReconstructX, uReconstructY, uReconstructZ, &
uReconstructZonal, uReconstructMeridional)
</font>
</pre>