<p><b>dwj07@fsu.edu</b> 2011-01-13 12:06:45 -0700 (Thu, 13 Jan 2011)</p><p><br>
        Fixing a bug with isopycnal coordinates in the ocean core and restart files.<br>
<br>
        Previously, Restart files were only read into the initial time slot, leaving<br>
        all other time levels initialized with zeros, to be updated during the simulations.<br>
        rho never gets updated in an isopycnal simulation, and thus would cause a NaN on the second<br>
        iteration. <br>
<br>
        Now, restart files are copied to all time levels, so all time levels have good data.<br>
<br>
        Also, added a feature to allow the use of config_stat_interval = 0, config_output_interval = 0 and<br>
        config_restart_interval = 0.<br>
<br>
        These options now supress the creation of statistic, output, and restart files, and don't cause mpas to error out anymore.<br>
</p><hr noshade><pre><font color="gray">Modified: trunk/mpas/src/core_ocean/module_mpas_core.F
===================================================================
--- trunk/mpas/src/core_ocean/module_mpas_core.F        2011-01-13 19:05:13 UTC (rev 683)
+++ trunk/mpas/src/core_ocean/module_mpas_core.F        2011-01-13 19:06:45 UTC (rev 684)
@@ -24,7 +24,6 @@
       type (block_type), pointer :: block
       type (dm_info) :: dminfo
 
-
       if (.not. config_do_restart) call setup_sw_test_case(domain)
 
       if (config_vert_grid_type.eq.'isopycnal') then
@@ -77,6 +76,7 @@
       type (block_type), intent(inout) :: block
       type (mesh_type), intent(inout) :: mesh
       real (kind=RKIND), intent(in) :: dt
+      integer :: i
    
    
       call compute_solve_diagnostics(dt, block % state % time_levs(1) % state, mesh)
@@ -85,7 +85,14 @@
       call init_reconstruct(mesh)
       call reconstruct(block % state % time_levs(1) % state, mesh)
    
-      if (.not. config_do_restart) block % state % time_levs(1) % state % xtime % scalar = 0.0
+      if (.not. config_do_restart) then 
+          block % state % time_levs(1) % state % xtime % scalar = 0.0
+      else
+          do i=2,nTimeLevs
+             call copy_state(block % state % time_levs(i) % state, &amp;
+                             block % state % time_levs(1) % state)
+          end do
+      endif
    
    end subroutine mpas_init_block
    
@@ -123,14 +130,18 @@
          ! Move time level 2 fields back into time level 1 for next time step
          call shift_time_levels_state(domain % blocklist % state)
    
-         if (mod(itimestep, config_output_interval) == 0) then
-            call write_output_frame(output_obj, output_frame, domain)
-         end if
-         if (mod(itimestep, config_restart_interval) == 0 .and. config_restart_interval &gt; 0) then
-            if (restart_frame == 1) call output_state_init(restart_obj, domain, &quot;RESTART&quot;)
-            call output_state_for_domain(restart_obj, domain, restart_frame)
-            restart_frame = restart_frame + 1
-         end if
+         if (config_output_interval &gt; 0) then
+             if (mod(itimestep, config_output_interval) == 0) then
+                 call write_output_frame(output_obj, output_frame, domain)
+             end if
+         endif
+         if (config_restart_interval &gt; 0) then
+             if (mod(itimestep, config_restart_interval) == 0) then
+                 if (restart_frame == 1) call output_state_init(restart_obj, domain, &quot;RESTART&quot;)
+                     call output_state_for_domain(restart_obj, domain, restart_frame)
+                     restart_frame = restart_frame + 1
+             end if
+         endif
       end do
 
    end subroutine mpas_core_run
@@ -208,18 +219,20 @@
    
       call timestep(domain, dt)
    
-      if (mod(itimestep, config_stats_interval) == 0) then
-         block_ptr =&gt; domain % blocklist
-         if(associated(block_ptr % next)) then
-            write(0,*) 'Error: computeGlobalDiagnostics assumes ',&amp;
-               'that there is only one block per processor.'
-         end if
+      if (config_stats_interval &gt; 0) then
+          if (mod(itimestep, config_stats_interval) == 0) then
+              block_ptr =&gt; domain % blocklist
+              if (associated(block_ptr % next)) then
+                  write(0,*) 'Error: computeGlobalDiagnostics assumes ',&amp;
+                            'that there is only one block per processor.'
+              end if
    
-         call timer_start(&quot;global diagnostics&quot;)
-         call computeGlobalDiagnostics(domain % dminfo, &amp;
-            block_ptr % state % time_levs(2) % state, block_ptr % mesh, &amp;
-            itimestep, dt)
-         call timer_stop(&quot;global diagnostics&quot;)
+          call timer_start(&quot;global diagnostics&quot;)
+          call computeGlobalDiagnostics(domain % dminfo, &amp;
+             block_ptr % state % time_levs(2) % state, block_ptr % mesh, &amp;
+             itimestep, dt)
+          call timer_stop(&quot;global diagnostics&quot;)
+          end if
       end if
    
    end subroutine mpas_timestep

</font>
</pre>