source: trunk/LMDZ.MARS/libf/phymars/nlthermeq.F @ 1747

Last change on this file since 1747 was 1266, checked in by aslmd, 11 years ago

LMDZ.MARS
IMPORTANT CHANGE

  • Remove all reference/use of nlayermx and dimphys.h
  • Made use of automatic arrays whenever arrays are needed with dimension nlayer
  • Remove lots of obsolete reference to dimensions.h
  • Converted iono.h and param_v4.h into corresponding modules

(with embedded subroutine to allocate arrays)
(no arrays allocated if thermosphere not used)

  • Deleted param.h and put contents into module param_v4_h
  • Adapted testphys1d, newstart, etc...
  • Made DATA arrays in param_read to be initialized by subroutine

fill_data_thermos in module param_v4_h

  • Optimized computations in paramfoto_compact (twice less dlog10 calculations)
  • Checked consistency before/after modification in debug mode
  • Checked performance is not impacted (same as before)
File size: 1.7 KB
Line 
1      subroutine nlthermeq(ngrid, nlayer, pplev, pplay)
2c
3c  Compute the number of layers nlaylte (stored in module yomlw_h)
4c  over which local thermodynamic equilibrium
5c  radiation scheme should be run to be sure of covering at least to a
6c  height greater than (pressure lower than) p=pminte, set in nlteparams.h.
7c  The maximum layer needed is found for the worst possible case.
8c  Stephen Lewis 6/2000
9c  Modified Y. Wanherdrick/ F. Forget 09/2000
10      use yomlw_h, only: nlaylte
11      implicit none
12#include "nlteparams.h"
13#include "callkeys.h"
14
15c
16c     Input:
17      integer ngrid, nlayer
18      real pplev(ngrid, nlayer+1)
19      real pplay(ngrid, nlayer)
20c
21c     Local:
22      integer igpmax, ismax
23      logical firstcall
24      data firstcall /.true./
25      save firstcall, igpmax
26c
27      if(firstcall) then
28c     Find the location of maximum surface pressure.
29c     Location won't vary much so only do it at the start;
30c     with no topography location would vary, but this is only
31c     needed for an estimate so any point would do in that case.
32         igpmax = ismax(ngrid, pplev, 1)
33         write(*, 10) ptrans
34         write(*, 20) zw
35         write(*, 30) pminte
36         firstcall = .false.
37      endif
38c
39      IF(callnlte) then
40c       Find first layer above pminte at this location
41        do nlaylte = nlayer, 1, -1
42           if (pplay(igpmax, nlaylte).gt.pminte)  go to 100
43        enddo
44      ELSE
45        nlaylte=nlayer
46      END IF
47  100 write(*,*) 'LTE rad. calculations up to layer ',  nlaylte
48c
49      return
50c
51   10 format(' nlthermeq: transition to NLTE centred at ',f6.2,'Pa')
52   20 format('               half-width (scale heights) ',f6.2)
53   30 format('          suggested LTE coverage at least ',f6.2,'Pa')
54      end
Note: See TracBrowser for help on using the repository browser.