<p><b>laura@ucar.edu</b> 2012-10-26 14:57:53 -0600 (Fri, 26 Oct 2012)</p><p>corrected the calculation of total pressure at cell vertices for interpolation of vorticity to fixed pressure levels<br>
</p><hr noshade><pre><font color="gray">Modified: branches/atmos_physics/src/core_nhyd_atmos/mpas_atm_interp_diagnostics.F
===================================================================
--- branches/atmos_physics/src/core_nhyd_atmos/mpas_atm_interp_diagnostics.F        2012-10-26 19:58:26 UTC (rev 2274)
+++ branches/atmos_physics/src/core_nhyd_atmos/mpas_atm_interp_diagnostics.F        2012-10-26 20:57:53 UTC (rev 2275)
@@ -35,7 +35,7 @@
real(kind=RKIND),dimension(:,:),pointer:: qvapor,theta_m,vorticity
real(kind=RKIND),dimension(:,:),pointer:: umeridional,uzonal,vvel
- real(kind=RKIND),dimension(:,:),allocatable:: pressure,pressure2,pressure_v,temperature
+ real(kind=RKIND),dimension(:,:),allocatable:: pressure,pressureCp1,pressure2,pressure_v,temperature
!local interpolated fields:
integer:: nIntP
@@ -48,7 +48,7 @@
write(0,*)
write(0,*) '--- enter subroutine interp_diagnostics:'
- nCells = mesh % nCellsSolve
+ nCells = mesh % nCells
nVertLevels = mesh % nVertLevels
!nVertLevelsP1 = mesh % nVertLevelsP1
nVertices = mesh % nVertices
@@ -71,17 +71,24 @@
umeridional => diag % uReconstructMeridional % array
uzonal => diag % uReconstructZonal % array
- if(.not.allocated(pressure) ) allocate(pressure(nVertLevels,nCells) )
- if(.not.allocated(pressure2) ) allocate(pressure2(nVertLevelsP1,nCells) )
- if(.not.allocated(pressure_v) ) allocate(pressure_v(nVertLevels,nVertices))
- if(.not.allocated(temperature)) allocate(temperature(nVertLevels,nCells) )
+ if(.not.allocated(pressure) ) allocate(pressure(nVertLevels,nCells) )
+ if(.not.allocated(pressureCp1) ) allocate(pressureCp1(nVertLevels,nCells+1) )
+ if(.not.allocated(pressure2) ) allocate(pressure2(nVertLevelsP1,nCells) )
+ if(.not.allocated(pressure_v) ) allocate(pressure_v(nVertLevels,nVertices) )
+ if(.not.allocated(temperature) ) allocate(temperature(nVertLevels,nCells) )
!calculation of total pressure at cell centers (at mass points):
do iCell = 1, nCells
do k = 1, nVertLevels
- pressure(k,iCell) = (pressure_p(k,iCell) + pressure_b(k,iCell)) / 100._RKIND
+ pressure(k,iCell) = (pressure_p(k,iCell) + pressure_b(k,iCell)) / 100._RKIND
+ pressureCp1(k,iCell) = pressure(k,iCell)
enddo
enddo
+ do iCell = nCells+1, nCells+1
+ do k = 1, nVertLevels
+ pressureCp1(k,iCell) = (pressure_p(k,iCell) + pressure_b(k,iCell)) / 100._RKIND
+ enddo
+ enddo
!calculation of total pressure at cell centers (at vertical velocity points):
k = nVertLevelsP1
@@ -117,8 +124,8 @@
do k = 1, nVertLevels
do iVertD = 1, vertexDegree
- pressure_v(k,iVert) = pressure_v(k,iVert) \
- + kiteAreasOnVertex(iVertD,iVert)*pressure(k,cellsOnVertex(iVertD,iVert))
+ pressure_v(k,iVert) = pressure_v(k,iVert) &
+ + kiteAreasOnVertex(iVertD,iVert)*pressureCp1(k,cellsOnVertex(iVertD,iVert))
enddo
pressure_v(k,iVert) = pressure_v(k,iVert) / areaTriangle(iVert)
enddo
@@ -270,6 +277,7 @@
if(allocated(field_interp)) deallocate(field_interp)
if(allocated(press_interp)) deallocate(press_interp)
if(allocated(pressure) ) deallocate(pressure )
+ if(allocated(pressureCp1) ) deallocate(pressureCp1 )
if(allocated(pressure2) ) deallocate(pressure2 )
if(allocated(pressure_v) ) deallocate(pressure_v )
if(allocated(temperature) ) deallocate(temperature )
</font>
</pre>