Ignore:
Timestamp:
Apr 3, 2014, 9:09:47 AM (11 years ago)
Author:
emillour
Message:

Generic model:
Major cleanup, in order to ease the use of LMDZ.GENERIC with (parallel) dynamics
in LMDZ.COMMON: (NB: this will break LMDZ.UNIVERSAL, which should be thrashed
in the near future)

  • Updated makegcm_* scripts (and makdim) and added the "-full" (to enforce full recomputation of the model) option
  • In dyn3d: converted control.h to module control_mod.F90 and converted iniadvtrac.F to module infotrac.F90
  • Added module mod_const_mpi.F90 in dyn3d (not used in serial mode)
  • Rearanged input/outputs routines everywhere to handle serial/MPI cases. physdem.F => phyredem.F90 , phyetat0.F => phyetat0.F90 ; all read/write routines for startfi files are gathered in module iostart.F90
  • added parallelism related routines init_phys_lmdz.F90, comgeomphy.F90, dimphy.F90, iniphysiq.F90, mod_grid_phy_lmdz.F90, mod_phys_lmdz_mpi_data.F90, mod_phys_lmdz_mpi_transfert.F90, mod_phys_lmdz_omp_data.F90, mod_phys_lmdz_omp_transfert.F90, mod_phys_lmdz_para.F90, mod_phys_lmdz_transfert_para.F90 in phymars and mod_const_mpi.F90 in dyn3d (for compliance with parallelism)
  • added created generic routines 'planetwide_maxval' and 'planetwide_minval', in module "planetwide_mod", that enable obtaining the max and min of a field over the whole planet. This should be further imroved with computation of means (possibly area weighed), etc.

EM

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.GENERIC/libf/phystd/soil_settings.F

    r787 r1216  
    1       subroutine soil_settings(nid,ngrid,nsoil,tsurf,tsoil)
    2 
    3       use comsoil_h
    4 
     1      subroutine soil_settings(nid,ngrid,nsoil,tsurf,tsoil,indextime)
     2
     3!      use netcdf
     4      use comsoil_h, only: layer, mlayer, inertiedat, volcapa
     5      use iostart, only: inquire_field_ndims, get_var, get_field,
     6     &                   inquire_field, inquire_dimension_length
    57      implicit none
    68
     
    911!
    1012!  Purpose: Read and/or initialise soil depths and properties
     13!
     14! Modifications: Aug.2010 EM : use NetCDF90 to load variables (enables using
     15!                      r4 or r8 restarts independently of having compiled
     16!                      the GCM in r4 or r8)
     17!                June 2013 TN : Possibility to read files with a time axis
     18!
    1119!
    1220!  This subroutine reads from a NetCDF file (opened by the caller)
     
    2634!======================================================================
    2735
    28 #include "dimensions.h"
    29 #include "dimphys.h"
    30 
    31 #include "netcdf.inc"
    3236!======================================================================
    3337!  arguments
    3438!  ---------
    3539!  inputs:
    36       integer nid       ! Input Netcdf file ID
    37       integer ngrid     ! # of horizontal grid points
    38       integer nsoil     ! # of soil layers
    39       integer sta       ! # at which reading starts
    40       real tsurf(ngrid) ! surface temperature
     40      integer,intent(in) :: nid ! Input Netcdf file ID
     41      integer,intent(in) :: ngrid       ! # of horizontal grid points
     42      integer,intent(in) :: nsoil       ! # of soil layers
     43      real,intent(in) :: tsurf(ngrid)   ! surface temperature
     44      integer,intent(in) :: indextime   ! position on time axis
    4145!  output:
    42       real tsoil(ngrid,nsoilmx) ! soil temperature
     46      real,intent(out) :: tsoil(ngrid,nsoil)    ! soil temperature
    4347
    4448!======================================================================
     
    5054      integer ndims     ! # of dimensions of read <inertiedat> data
    5155      integer ig,iloop  ! loop counters
     56     
     57      integer edges(3),corner(3) ! to read a specific time
    5258
    5359      logical :: olddepthdef=.false. ! flag
     
    7076      real,parameter :: default_volcapa=1.e6
    7177
    72       integer, dimension(2) :: sta2d
    73       integer, dimension(2) :: siz2d
    74 
    75 !======================================================================
    76 
    77       !! this is defined in comsoil_h
    78       IF (.not.ALLOCATED(layer))
    79      . ALLOCATE(layer(nsoilmx))
    80       IF (.not.ALLOCATED(mlayer))
    81      . ALLOCATE(mlayer(0:nsoilmx-1))
    82       IF (.not.ALLOCATED(inertiedat))
    83      . ALLOCATE(inertiedat(ngrid,nsoilmx))
    84 
    85       !! this is defined in dimphys.h
    86       sta = cursor
     78      logical :: found,ok
     79     
     80!======================================================================
    8781
    8882! 1. Depth coordinate
     
    9185! 1.1 Start by reading how many layers of soil there are
    9286
    93         ierr=NF_INQ_DIMID(nid,"subsurface_layers",dimid)
    94         if (ierr.ne.NF_NOERR) then
    95          write(*,*)'soil_settings: Problem reading <subsurface_layers>'
    96          call abort
    97         endif
    98 
    99         ierr=NF_INQ_DIMLEN(nid,dimid,dimlen)
    100         if (ierr.ne.NF_NOERR) then
    101          write(*,*)'soil_settings: Problem getting <subsurface_layers>',
    102      &             'length'
    103          call abort
    104         endif
     87        dimlen=inquire_dimension_length("subsurface_layers")
    10588
    10689        if (dimlen.ne.nsoil) then
     
    11295        endif
    11396
    114       sta2d = (/sta,1/)
    115       siz2d = (/ngrid,dimlen/)
    116 
    117 
    11897! 1.2 Find out the # of dimensions <inertiedat> was defined as using
    119 !     (in ye old days, thermal inertia was only given at the "surface")
    120       ! Look for thermal inertia data
    121       ierr=NF_INQ_VARID(nid, "inertiedat", nvarid)
    122       if (ierr.NE.NF_NOERR) then
    123          write(*,*)'soil_settings: Field <inertiedat> not found!'
    124          call abort
    125       endif
    126 
    127       ! Read the # of dimensions <inertidat> was defined as using
    128       ierr=NF_INQ_VARNDIMS(nid,nvarid,ndims)
    129       ! if (ndims.eq.1) then we have the "old 2D-surface" format
     98
     99      ndims=inquire_field_ndims("inertiedat")
    130100
    131101! 1.3 Read depths values or set olddepthdef flag and values
     
    147117        enddo
    148118      else ! Look for depth
    149         ierr=NF_INQ_VARID(nid,"soildepth",nvarid)
    150         if (ierr.ne.NF_NOERR) then
    151           write(*,*)'soil_settings: Field <soildepth> not found!'
    152           call abort
    153         endif
    154119        ! read <depth> coordinate
    155120        if (interpol) then !put values in oldmlayer
    156           if (.not.allocated(oldmlayer)) then
    157              allocate(oldmlayer(dimlen),stat=ierr)
    158           endif
    159 #ifdef NC_DOUBLE
    160            ierr = NF_GET_VAR_DOUBLE(nid,nvarid,oldmlayer)
    161 #else
    162            ierr = NF_GET_VAR_REAL(nid,nvarid,oldmlayer)
    163 #endif
    164           if (ierr.ne.NF_NOERR) then
    165            write(*,*)'soil_settings: Problem while reading <soildepth>'
    166            call abort
     121          call get_var("soildepth",oldmlayer,found)
     122          if (.not.found) then
     123            write(*,*)'soil_settings: Problem while reading <soildepth>'
    167124          endif
    168125        else ! put values in mlayer
    169 #ifdef NC_DOUBLE
    170            ierr = NF_GET_VAR_DOUBLE(nid,nvarid,mlayer)
    171 #else
    172            ierr = NF_GET_VAR_REAL(nid,nvarid,mlayer)
    173 #endif
    174           if (ierr.ne.NF_NOERR) then
    175            write(*,*)'soil_settings: Problem while reading <soildepth>'
    176            call abort
     126          call get_var("soildepth",mlayer,found)
     127          if (.not.found) then
     128            write(*,*)'soil_settings: Problem while reading <soildepth>'
    177129          endif
    178130        endif !of if (interpol)
     
    183135      ! default mlayer distribution, following a power law:
    184136      !  mlayer(k)=lay1*alpha**(k-1/2)
    185         lay1=2.e-4  !mars case (nsoilmx=18)
    186 !        lay1=3.e-2  !earth case (nsoilmx=13)
     137        lay1=2.e-4
    187138        alpha=2
    188139        do iloop=0,nsoil-1
     
    200151      enddo
    201152
    202 ! 2. Volumetric heat capacity (note: it is declared in comsoil.h)
     153! 2. Volumetric heat capacity (note: it is declared in comsoil_h)
    203154! ---------------------------
    204155! "volcapa" is (so far) 0D and written in "controle" table of startfi file
     
    214165        volcapa=default_volcapa
    215166      endif
    216 ! Look for it
    217 !      ierr=NF_INQ_VARID(nid,"volcapa",nvarid)
    218 !      if (ierr.NE.NF_NOERR) then
    219 !         write(*,*)'soil_settings: Field <volcapa> not found!'
    220 !         write(*,*)'Initializing Volumetric heat capacity to ',
    221 !     &             default_volcapa
    222 !         volcapa=default_volcapa
    223 !      else
    224 !#ifdef NC_DOUBLE
    225 !       ierr = NF_GET_VAR_DOUBLE(nid,nvarid,volcapa)
    226 !#else
    227 !       ierr = NF_GET_VAR_REAL(nid,nvarid,volcapa)
    228 !#endif
    229 !        if (ierr.ne.NF_NOERR) then
    230 !         write(*,*)'soil_settings: Problem while reading <volcapa>'
    231 !         call abort
    232 !       endif
    233 !      endif
    234 
    235 ! 3. Thermal inertia (note: it is declared in comsoil.h)
     167
     168! 3. Thermal inertia (note: it is declared in comsoil_h)
    236169! ------------------
    237170
    238171! 3.1 Look (again) for thermal inertia data (to reset nvarid)
    239       ierr=NF_INQ_VARID(nid, "inertiedat", nvarid)
    240       if (ierr.NE.NF_NOERR) then
    241          write(*,*)'soil_settings: Field <inertiedat> not found!'
    242          call abort
    243       endif
    244172
    245173! 3.2 Knowing the # of dimensions <inertidat> was defined as using,
     
    251179       ! Read Surface thermal inertia
    252180       allocate(surfinertia(ngrid))
    253 #ifdef NC_DOUBLE
    254        ierr = NF_GET_VARA_DOUBLE(nid, nvarid, sta, ngrid, surfinertia)
    255 #else
    256        ierr = NF_GET_VARA_REAL(nid, nvarid, sta, ngrid, surfinertia)
    257 #endif
    258         if (ierr.NE.NF_NOERR) then
    259          write(*,*)'soil_settings: Problem while reading <inertiedat>'
     181       call get_field("inertiedat",surfinertia,found)
     182       if (.not.found) then
     183         write(*,*) "soil_settings: Failed loading <inertiedat>"
    260184         call abort
    261         endif
     185       endif
    262186       
    263187       write(*,*)' => Building soil thermal inertia (using reference sur
     
    277201            stop
    278202           endif
    279          endif
    280 #ifdef NC_DOUBLE
    281         ierr = NF_GET_VARA_DOUBLE(nid,nvarid,sta2d,siz2d,oldinertiedat)
    282 #else
    283         ierr = NF_GET_VARA_REAL(nid,nvarid,sta2d,siz2d,oldinertiedat)
    284 #endif
    285         if (ierr.NE.NF_NOERR) then
    286          write(*,*)'soil_settings: Problem while reading <inertiedat>'
    287          call abort
     203         endif ! of if (.not.allocated(oldinertiedat))
     204        call get_field("inertiedat",oldinertiedat,found)
     205        if (.not.found) then
     206          write(*,*) "soil_settings: Failed loading <inertiedat>"
     207          call abort
    288208        endif
    289209       else ! put values in therm_i
    290 #ifdef NC_DOUBLE
    291         ierr = NF_GET_VARA_DOUBLE(nid,nvarid,sta2d,siz2d,inertiedat)
    292 #else
    293         ierr = NF_GET_VARA_REAL(nid,nvarid,sta2d,siz2d,inertiedat)
    294 #endif
    295         if (ierr.NE.NF_NOERR) then
    296          write(*,*)'soil_settings: Problem while reading <inertiedat>'
    297          call abort
    298         endif
     210         call get_field("inertiedat",inertiedat,found)
     211         if (.not.found) then
     212           write(*,*) "soil_settings: Failed loading <inertiedat>"
     213           call abort
     214         endif
     215!        endif
    299216       endif ! of if (interpol)
    300217      endif ! of if (ndims.eq.1)
     
    303220! -------------------------
    304221
    305       ierr=NF_INQ_VARID(nid,"tsoil",nvarid)
    306       if (ierr.ne.NF_NOERR) then
     222!      ierr=nf90_inq_varid(nid,"tsoil",nvarid)
     223      ok=inquire_field("tsoil")
     224!      if (ierr.ne.nf90_noerr) then
     225      if (.not.ok) then
    307226        write(*,*)'soil_settings: Field <tsoil> not found!'
    308227        write(*,*)' => Building <tsoil> from surface values <tsurf>'
     
    320239           endif
    321240         endif
    322 #ifdef NC_DOUBLE
    323         ierr=NF_GET_VARA_DOUBLE(nid,nvarid,sta2d,siz2d,oldtsoil)
    324 #else
    325         ierr=NF_GET_VARA_REAL(nid,nvarid,sta2d,siz2d,oldtsoil)
    326 #endif
    327         if (ierr.ne.NF_NOERR) then
    328          write(*,*)'soil_settings: Problem while reading <tsoil>'
    329          call abort
    330         endif
     241         call get_field("tsoil",oldtsoil,found)
     242         if (.not.found) then
     243           write(*,*) "soil_settings: Failed loading <tsoil>"
     244           call abort
     245         endif
    331246       else ! put values in tsoil
    332 #ifdef NC_DOUBLE
    333         ierr=NF_GET_VARA_DOUBLE(nid,nvarid,sta2d,siz2d,tsoil)
    334 #else
    335         ierr=NF_GET_VARA_REAL(nid,nvarid,sta2d,siz2d,tsoil)
    336 #endif
    337         if (ierr.ne.NF_NOERR) then
    338          write(*,*)'soil_settings: Problem while reading <tsoil>'
    339          call abort
    340         endif
     247         call get_field("tsoil",tsoil,found,timeindex=indextime)
     248         if (.not.found) then
     249           write(*,*) "soil_settings: Failed loading <tsoil>"
     250           call abort
     251         endif
    341252       endif ! of if (interpol)
    342       endif
     253      endif! of if (.not.ok)
    343254
    344255! 5. If necessary, interpolate soil temperatures and thermal inertias
Note: See TracChangeset for help on using the changeset viewer.