updated top-level README and version_decl for V4.5 (#1847)
[WRF.git] / share / module_HLaw.F
blob5363d7ff2260528d10a1be0adb1098106ae629d5
2    module module_HLawConst
4    implicit none
6    private
8    integer :: nHLC
10    type HLCnst_type
11      real              :: mw
12      real              :: hcnst(6)
13      character(len=64) :: name
14    end type HLCnst_type
16    type(HLCnst_type), allocatable :: HLC(:)
18    public :: init_HLawConst
19    public :: HLCnst_type
20    public :: HLC, nHLC
22    CONTAINS
24    subroutine init_HLawConst( dm )
26    integer, intent(in) :: dm                  ! domain index
28    integer :: m, n, unitno
29    integer :: nNull
30    integer :: astat, istat
31    integer :: LawType
32    real    :: inMw
33    real    :: inHeff(6)
34    character(len=64)  :: inName
35    character(len=256) :: inlin
36    character(len=256) :: emsg
38    integer, external :: get_unused_unit
40 top_level_domain: &
41    if( dm == 1 .and. .not. allocated(HLC) ) then
42      unitno = get_unused_unit()
43      if( unitno <= 0 ) then
44        call wrf_error_fatal( 'init_HLConst: Failed to get Fortran I/O unit number' )
45      endif
46      open(unit=unitno,file='HLC.TBL',status='OLD',iostat=istat)
47      if( istat /= 0 ) then
48        write(emsg,'(''init_HLConst: Failed to open HLC.TBL; error = '',i8)') istat
49        call wrf_error_fatal( trim(emsg) )
50      endif
52      nHLC = 0
53      do
54        read(unitno,*,iostat=istat) inlin
55        if( istat /= 0 ) then
56          exit
57        else
58          nHLC = nHLC + 1
59        endif
60      end do
62      write(emsg,'(''Read '',i4,'' Henrys Law entries'')') nHLC
63      call wrf_debug( 0,trim(emsg) )
65      if( nHLC > 0 ) then
66        allocate( HLC(nHLC),stat=astat )
67        if( astat /= 0 ) then 
68          write(emsg,'(''init_HLawConst: Failed to allocate HLC; error = '',i8)') astat
69          call wrf_error_fatal( trim(emsg) )
70        endif
71        rewind(unit=unitno)
72        nNull = 0 ; m = 0
73        do n = 1,nHLC
74          read(unitno,*,iostat=istat) inName,LawType,inMw,inHeff
75          if( istat /= 0 ) then
76            write(emsg,'(''init_HLawConst: Failed to read line '',i3,''; error = '',i8)') n,istat
77            call wrf_error_fatal( trim(emsg) )
78          else
79            if( all( inHeff == 0. ) ) then
80              nNull = nNull + 1
81              cycle
82            else
83              m = m + 1
84              HLC(m)%name = inName ; HLC(m)%mw = inMw ; HLC(m)%hcnst(:) = inHeff(:)
85            endif
86          endif
87        end do
89        write(emsg,'(''There are '',i4,'' Henrys Law null entries'')') nNull
90        call wrf_debug( 0,trim(emsg) )
91        nHLC = nHLC - nNull
92        call wrf_debug( 0, ' ' )
93        call wrf_debug( 0, 'HLaw table ' )
94        do n = 1,nHLC
95          write(emsg,'(''('',i3.3,'')'',a16,1pg15.7,3x,6g15.7)') &
96             n,trim(HLC(n)%name),HLC(n)%mw,HLC(n)%hcnst(:)
97          call wrf_debug( 0, trim(emsg) )
98        end do
99      endif
101      close(unit=unitno)
103    endif top_level_domain
105    end subroutine init_HLawConst
107    end module module_HLawConst