Opened 5 years ago

Closed 2 years ago

#140 closed Defect (fixed)

vertransform.f90 occasionally fails to initialize

Reported by: adingwell Owned by: pesei
Priority: major Milestone: FLEXPART 10
Component: FP other Version: FLEXPART 9.2
Keywords: Cc:

Description (last modified by pesei)

When running with smaller domains verttransform might fail to initialize ixm and jym.

At initialization the following block is run in verttransform.f90:

  ! Search for a point with high surface pressure (i.e. not above significant topography)
  ! Then, use this point to construct a reference z profile, to be used at all times
  !**************************************************************************************

    do jy=0,nymin1
      do ix=0,nxmin1
        if (ps(ix,jy,1,n).gt.100000.) then
          ixm=ix
          jym=jy
          goto 3
        endif
      end do
    end do
3   continue

In order to work, there has to be a grid point within the domain with a surface pressure above 1000 hPa. Such a point might not be present in smaller domains, where there is only a small number of grid points over sea surface.

It would be good if this block at least produces an error if no reference point is found.

Change History (5)

comment:1 Changed 5 years ago by pesei

  • Owner set to pesei
  • Status changed from new to accepted

I suggest that we take the lowest grid point found if there is none >1000 hPa.
I'll take the ticket for v9.2.

comment:2 follow-up: Changed 2 years ago by adingwell

I just noticed that this bug found its way into 10.2.

I've added the following lines between the loop and the continue statement in my local version:

ij_min  = MINLOC(ps(:,:,1,n))  ! ij-index of minimum pressure
ixm = ij_min(1)
jym = ij_min(2)

but I haven't tested it for a problematic domain yet (as my problems this time were related to #197)

Last edited 2 years ago by pesei (previous) (diff)

comment:3 Changed 2 years ago by pesei

  • Description modified (diff)
  • Milestone changed from FLEXPART 9.2 to FLEXPART 10

comment:4 in reply to: ↑ 2 Changed 2 years ago by pesei

Replying to adingwell:

I just noticed that this bug found its way into 10.2.

I've added the following lines between the loop and the continue statement in my local version:

ij_min  = MINLOC(ps(:,:,1,n))  ! ij-index of minimum pressure

don't use min it is the opposite of what we do in the code! We have to use maxloc to find the point the ca lowest elevation. Also, we have to be aware that this will cause results to differ from precvous versions, so maybe either offer option for old one or do this only if the old condition fails.

comment:5 Changed 2 years ago by pesei

  • Resolution set to fixed
  • Status changed from accepted to closed

This is now fixed in dev, see
verttransform_ecmwf: changeset:1a8fbee08b4dd8d4aab8d4f78b6bd0bb37e21f13/flexpart.git
verttransform_gfs: changeset:f251e576e36228966cea4e1a90c93e625f5cd2d9/flexpart.git

Currently, there is a parameter in subroutine header which can be set to 0 for the old selection of the reference grid cell, or 1 for the new method. This might be moved to par_mod.f90 later. The method selects the grid point whose sfc pressure is closest to the mean over the domain (skipping regions with pressure lower than first-guess mean plus 1 sigma).

For a global GFS data set, this is a significant difference. The old method takes a very cold Antarctic grid cell and has about 10% lower heights.

Note: See TracTickets for help on using tickets.
hosted by ZAMG