Ignore:
Timestamp:
Mar 14, 2023, 10:07:33 AM (22 months ago)
Author:
romain.vande
Message:

Mars PCM:
Adapt start2archive.F to the subslope parametrisation.
Small correction for some dimensions of variables.
RV

Location:
trunk/LMDZ.MARS/libf/phymars
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.MARS/libf/phymars/iostart.F90

    r2907 r2913  
    1515    INTEGER,SAVE :: idim7 ! "Time" dimension
    1616    INTEGER,SAVE :: idim8 ! "nslope" dimension
    17     INTEGER,SAVE :: idim9 ! "inter slope" dimension
     17    INTEGER,SAVE :: idim9 ! "nslope_plus_1" dimension
    1818    INTEGER,SAVE :: timeindex ! current time index (for time-dependent fields)
    1919    INTEGER,PARAMETER :: length=100 ! size of tab_cntrl array
    … …  
    564564      ENDIF
    565565
    566       ierr=NF90_DEF_DIM(nid_restart,"inter slope",nslope+1,idim9)
     566      ierr=NF90_DEF_DIM(nid_restart,"inter_slope",nslope+1,idim9)
    567567      IF (ierr/=NF90_NOERR) THEN
    568568        write(*,*)'phyredem: problem defining inter slope dimension'
    … …  
    660660 
    661661  REAL                           :: field_glo(klon_glo,field_size)
     662  REAL                           :: field_glo_reshape(klon_glo,nsoilmx,nslope,timeindex)
    662663  INTEGER                        :: ierr
    663664  INTEGER                        :: nvarid
    … …  
    804805        endif ! of if (.not.present(timeindex))
    805806
    806       ELSE IF (field_size==nsoilmx) THEN
     807      ELSE IF (field_size==nsoilmx .or. field_size==nsoilmx*nslope) THEN
    807808        ! input is a 2D "subsurface field" array
    808809        if (.not.present(time)) then ! for a time-independent field
    … …  
    828829            ierr=NF90_REDEF(nid_restart)
    829830#ifdef NC_DOUBLE
     831            if(field_name.eq. "tsoil") then
     832            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
     833                              (/idim2,idim3,idim8,idim7/),nvarid)
     834            else
    830835            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,&
    831836                              (/idim2,idim3,idim7/),nvarid)
    832 #else
     837            endif
     838#else
     839            if(field_name.eq. "tsoil") then
     840            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
     841                              (/idim2,idim3,idim8,idim7/),nvarid)
     842            else
    833843            ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,&
    834844                              (/idim2,idim3,idim7/),nvarid)
     845            endif
    835846#endif
    836847           if (ierr.ne.NF90_NOERR) then
    … …  
    842853          endif
    843854          ! Write the variable
     855          if(field_name.eq. "tsoil") then
     856
     857          field_glo_reshape=RESHAPE(field_glo, (/klon_glo, nsoilmx, nslope, timeindex/))
     858
     859          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo_reshape,&
     860                            start=(/1,1,1,timeindex/))
     861          else
    844862          ierr=NF90_PUT_VAR(nid_restart,nvarid,field_glo,&
    845863                            start=(/1,1,timeindex/))
    846 
     864          endif
    847865        endif ! of if (.not.present(time))
    848866

  • trunk/LMDZ.MARS/libf/phymars/phyetat0_mod.F90

    r2909 r2913  
    6060
    6161!  outputs:
    62   real,intent(out) :: tsurf(ngrid) ! surface temperature
    63   real,intent(out) :: tsoil(ngrid,nsoil) ! soil temperature
    64   real,intent(out) :: albedo(ngrid,2) ! surface albedo
    65   real,intent(out) :: emis(ngrid) ! surface emissivity
     62  real,intent(out) :: tsurf(ngrid,nslope) ! surface temperature
     63  real,intent(out) :: tsoil(ngrid,nsoil,nslope) ! soil temperature
     64  real,intent(out) :: albedo(ngrid,2,nslope) ! surface albedo
     65  real,intent(out) :: emis(ngrid,nslope) ! surface emissivity
    6666  real,intent(out) :: q2(ngrid,nlay+1) !
    67   real,intent(out) :: qsurf(ngrid,nq) ! tracers on surface
     67  real,intent(out) :: qsurf(ngrid,nq,nslope) ! tracers on surface
    6868  real,intent(out) :: tauscaling(ngrid) ! dust conversion factor
    6969  real,intent(out) :: totcloudfrac(ngrid) ! total cloud fraction
    7070  real,intent(out) :: wstar(ngrid) ! Max vertical velocity in thermals (m/s)
    71   real,intent(out) :: watercap(ngrid) ! h2o_ice_cover
     71  real,intent(out) :: watercap(ngrid,nslope) ! h2o_ice_cover
    7272  real, intent(out) :: def_slope(nslope+1) !boundaries for bining of the slopes
    7373  real, intent(out) :: def_slope_mean(nslope)
    … …  
    143143      write(*,*)'We take default def_slope and subslope dist'
    144144      allocate(default_def_slope(nslope+1))
    145       default_def_slope(1) = 0.
    146       default_def_slope(2) = 0.
     145      default_def_slope(1) = -90.
     146      default_def_slope(2) = 90.
    147147      subslope_dist(:,:)=1.
     148      def_slope(:)=default_def_slope(:)
    148149     else
    149150       write(*,*)'Variable def_slope is not present in the start and'
    … …  
    168169      default_def_slope(2) = 0.
    169170      subslope_dist(:,:)=1.
     171      def_slope(:)=default_def_slope(:)
    170172    else
    171173     write(*,*)'Without starfi, lets first run with nslope=1'
    … …  
    516518   endif
    517519else
    518   tsurf(:)=0. ! will be updated afterwards in physiq !
     520  tsurf(:,:)=0. ! will be updated afterwards in physiq !
    519521endif ! of if (startphy_file)
    520522write(*,*) "phyetat0: Surface temperature <tsurf> range:", &
    … …  
    523525! Surface albedo
    524526if (startphy_file) then
    525    call get_field("albedo",albedo(:,1),found,indextime)
     527   call get_field("albedo",albedo(:,1,:),found,indextime)
    526528   if (.not.found) then
    527529     write(*,*) "phyetat0: Failed loading <albedo>"
    528      albedo(:,1)=albedodat(:)
    529    endif
    530 else
    531    albedo(:,1)=albedodat(:)
     530     do islope=1,nslope
     531       albedo(:,1,islope)=albedodat(:)
     532     enddo
     533   endif
     534else
     535   do islope=1,nslope
     536     albedo(:,1,islope)=albedodat(:)
     537   enddo
    532538endif ! of if (startphy_file)
    533539write(*,*) "phyetat0: Surface albedo <albedo> range:", &
    534             minval(albedo(:,1)), maxval(albedo(:,1))
    535 albedo(:,2)=albedo(:,1)
     540            minval(albedo(:,1,:)), maxval(albedo(:,1,:))
     541albedo(:,2,:)=albedo(:,1,:)
    536542
    537543! Surface emissivity
    … …  
    547553  call getin_p("surfemis_without_startfi",surfemis)
    548554  print*,"surfemis_without_startfi",surfemis
    549   emis(:)=surfemis
     555  emis(:,:)=surfemis
    550556endif ! of if (startphy_file)
    551557write(*,*) "phyetat0: Surface emissivity <emis> range:", &
    … …  
    633639        ! We first check if co2ice exist in the startfi file (old way)
    634640        ! CO2 ice cover
    635         call get_field("co2ice",qsurf(:,iq),found,indextime)
     641        call get_field("co2ice",qsurf(:,iq,:),found,indextime)
    636642        ! If not, we load the corresponding tracer in qsurf
    637643        if (.not.found) then
    638           call get_field(txt,qsurf(:,iq),found,indextime)
     644          call get_field(txt,qsurf(:,iq,:),found,indextime)
    639645          if (.not.found) then
    640646            call abort_physic(modname, &
    … …  
    643649        endif
    644650      else ! (not the co2 tracer)
    645         call get_field(txt,qsurf(:,iq),found,indextime)
     651        call get_field(txt,qsurf(:,iq,:),found,indextime)
    646652        if (.not.found) then
    647653          write(*,*) "phyetat0: Failed loading <",trim(txt),">"
    648654          write(*,*) "         ",trim(txt)," is set to zero"
    649           qsurf(:,iq)=0.
     655          qsurf(:,iq,:)=0.
    650656        endif
    651657      endif !endif co2
    652658    else !(not startphy_file)
    653       qsurf(:,iq)=0.
     659      qsurf(:,iq,:)=0.
    654660    endif ! of if (startphy_file)
    655661    write(*,*) "phyetat0: Surface tracer <",trim(txt),"> range:", &
    656                  minval(qsurf(:,iq)), maxval(qsurf(:,iq))
     662                 minval(qsurf(:,iq,:)), maxval(qsurf(:,iq,:))
    657663  enddo ! of do iq=1,nq
    658664
    … …  
    661667      ! CO2 ice cover
    662668      if (startphy_file) then
    663         call get_field("co2ice",qsurf(:,iq),found,indextime)
     669        call get_field("co2ice",qsurf(:,iq,:),found,indextime)
    664670      ! If not, we load the corresponding tracer in qsurf
    665671        if (.not.found) then
    666           call get_field(txt,qsurf(:,iq),found,indextime)
     672          call get_field(txt,qsurf(:,iq,:),found,indextime)
    667673          if (.not.found) then
    668674            call abort_physic(modname, &
    … …  
    672678      else
    673679       ! If we run without startfile, co2ice is set to 0
    674         qsurf(:,iq)=0.
     680        qsurf(:,iq,:)=0.
    675681      endif !if (startphy_file)
    676682        write(*,*) "phyetat0: CO2 ice cover <co2ice> range:", &
    677             minval(qsurf(:,iq)), maxval(qsurf(:,iq))
     683            minval(qsurf(:,iq,:)), maxval(qsurf(:,iq,:))
    678684    endif
    679685
    … …  
    685691     write(*,*) "phyetat0: Failed loading <watercap> : ", &
    686692                          "<watercap> is set to zero"
    687      watercap(:)=0
     693     watercap(:,:)=0
    688694     write(*,*) 'Now transfer negative surface water ice to', &
    689695                ' watercap !'
    … …  
    693699       if (txt.eq."h2o_ice") then
    694700         do ig=1,ngrid
    695           if (qsurf(ig,iq).lt.0.0) then
    696              watercap(ig) = qsurf(ig,iq)
    697              qsurf(ig,iq) = 0.0
    698           end if
     701          do islope=1,nslope
     702            if (qsurf(ig,iq,islope).lt.0.0) then
     703              watercap(ig,islope) = qsurf(ig,iq,islope)
     704              qsurf(ig,iq,islope) = 0.0
     705            end if
     706          enddo
    699707         end do
    700708       endif
    … …  
    702710       if (txt.eq."hdo_ice") then
    703711         do ig=1,ngrid
    704           if (qsurf(ig,iq).lt.0.0) then
    705              qsurf(ig,iq) = 0.0
    706           end if
     712          do islope=1,nslope
     713            if (qsurf(ig,iq,islope).lt.0.0) then
     714              qsurf(ig,iq,islope) = 0.0
     715            end if
     716          enddo
    707717         end do
    708718       endif
    … …  
    712722   endif ! of if (.not.found)
    713723else
    714    watercap(:)=0
     724   watercap(:,:)=0
    715725endif ! of if (startphy_file)
    716726write(*,*) "phyetat0: Surface water ice <watercap> range:", &

  • trunk/LMDZ.MARS/libf/phymars/phyredem.F90

    r2896 r2913  
    171171  use dust_param_mod, only: dustscaling_mode
    172172  use comsoil_h,only: flux_geo
     173  use comslope_mod, ONLY: nslope
    173174  implicit none
    174175 
    … …  
    182183  real,intent(in) :: phystep
    183184  real,intent(in) :: time
    184   real,intent(in) :: tsurf(ngrid)
    185   real,intent(in) :: tsoil(ngrid,nsoil)
    186   real,intent(in) :: albedo(ngrid,2)
    187   real,intent(in) :: emis(ngrid)
     185  real,intent(in) :: tsurf(ngrid,nslope)
     186  real,intent(in) :: tsoil(ngrid,nsoil,nslope)
     187  real,intent(in) :: albedo(ngrid,2,nslope)
     188  real,intent(in) :: emis(ngrid,nslope)
    188189  real,intent(in) :: q2(ngrid,nlay+1)
    189   real,intent(in) :: qsurf(ngrid,nq)
     190  real,intent(in) :: qsurf(ngrid,nq,nslope)
    190191  real,intent(in) :: tauscaling(ngrid)
    191192  real,intent(in) :: totcloudfrac(ngrid)
    192193  real,intent(in) :: wstar(ngrid)
    193   real,intent(in) :: watercap(ngrid)
     194  real,intent(in) :: watercap(ngrid,nslope)
    194195 
    195196  integer :: iq
    … …  
    219220 
    220221  ! Albedo of the surface
    221   call put_field("albedo","Surface albedo",albedo(:,1),time)
     222  call put_field("albedo","Surface albedo",albedo(:,1,:),time)
    222223 
    223224  ! Emissivity of the surface
    … …  
    312313      end if     
    313314
    314       call put_field(trim(txt),"tracer on surface",qsurf(:,iq),time)
     315      call put_field(trim(txt),"tracer on surface",qsurf(:,iq,:),time)
    315316    enddo
    316317  endif

  • trunk/LMDZ.MARS/libf/phymars/soil_settings.F

    r2887 r2913  
    55      use iostart, only: inquire_field_ndims, get_var, get_field,
    66     &                   inquire_field, inquire_dimension_length
     7      USE comslope_mod, ONLY: nslope
    78      implicit none
    89
    … …  
    4142      integer,intent(in) :: ngrid       ! # of horizontal grid points
    4243      integer,intent(in) :: nsoil       ! # of soil layers
    43       real,intent(in) :: tsurf(ngrid)   ! surface temperature
     44      real,intent(in) :: tsurf(ngrid,nslope)   ! surface temperature
    4445      integer,intent(in) :: indextime   ! position on time axis
    4546!  output:
    46       real,intent(out) :: tsoil(ngrid,nsoil)    ! soil temperature
     47      real,intent(out) :: tsoil(ngrid,nsoil,nslope)     ! soil temperature
    4748
    4849!======================================================================
    … …  
    5455      integer ndims     ! # of dimensions of read <inertiedat> data
    5556      integer ig,iloop  ! loop counters
     57      integer islope
    5658     
    5759      integer edges(3),corner(3) ! to read a specific time
    … …  
    6466      real,dimension(:),allocatable :: oldmlayer
    6567      real,dimension(:,:),allocatable :: oldinertiedat
    66       real,dimension(:,:),allocatable :: oldtsoil
     68      real,dimension(:,:,:),allocatable :: oldtsoil
    6769     
    6870      ! for interpolation
    … …  
    317319        write(*,*)' => Building <tsoil> from surface values <tsurf>'
    318320        do iloop=1,nsoil
    319           tsoil(:,iloop)=tsurf(:)
     321          do islope=1,nslope
     322            tsoil(:,iloop,islope)=tsurf(:,islope)
     323          enddo
    320324        enddo
    321325      else ! <tsoil> found
    322326       if (interpol) then ! put values in oldtsoil
    323327         if (.not.allocated(oldtsoil)) then
    324            allocate(oldtsoil(ngrid,dimlen),stat=ierr)
     328           allocate(oldtsoil(ngrid,dimlen,nslope),stat=ierr)
    325329           if (ierr.ne.0) then
    326330             write(*,*) 'soil_settings: failed allocation of ',
    … …  
    376380        do ig=1,ngrid
    377381          ! copy values
    378           oldval(1)=tsurf(ig)
    379           oldval(2:dimlen+1)=oldtsoil(ig,1:dimlen)
     382          do islope=1,nslope
     383          oldval(1)=tsurf(ig,islope)
     384          oldval(2:dimlen+1)=oldtsoil(ig,1:dimlen,islope)
    380385          ! build vertical coordinate
    381386          oldgrid(1)=0. ! ground
    … …  
    385390          call interp_line(oldgrid,oldval,dimlen+1,mlayer,newval,nsoil)
    386391          ! copy result in tsoil
    387           tsoil(ig,:)=newval(:)
     392          tsoil(ig,:,islope)=newval(:)
     393          enddo
    388394        enddo
    389395
    … …  
    419425        oldgrid(2:dimlen+1)=oldmlayer(1:dimlen)
    420426        do ig=1,ngrid
     427          do islope=1,nslope
    421428          ! data
    422           oldval(1)=tsurf(ig)
    423           oldval(2:dimlen+1)=oldtsoil(ig,1:dimlen)
     429          oldval(1)=tsurf(ig,islope)
     430          oldval(2:dimlen+1)=oldtsoil(ig,1:dimlen,islope)
    424431          ! interpolate
    425432          call interp_line(oldgrid,oldval,dimlen+1,mlayer,newval,nsoil)
    426433          ! copy result in inertiedat
    427           tsoil(ig,:)=newval(:)
     434          tsoil(ig,:,islope)=newval(:)
     435          enddo
    428436        enddo
    429437       
Note: See TracChangeset for help on using the changeset viewer.