!=================================================
! 3D cumulus model - kaminari
!
! FileName : base.f90
! Note     : $B4pK\>l7W;;%W%m%0%i%`(B
!            $BC10L$O(B MKS $B7O(B
!            $BL5<!8505NO(B = $B%(%/%9%J!<4X?t(B
!            $B$^$:(B 2$B<!85%b%G%k$r:n$k(B (x-z $B7O(B)
! Author   : TAKAHASHI Koko
! Date     : 2003/10/30
!            2003/10/29
!            2003/09/01 $B?75,:n@.(B
! 
!=================================================
program base

  !-------------------------------------------------
  ! $BJQ?t$J$I$N@_Dj(B
  !-------------------------------------------------

  use gtool_history
  implicit none

  ! $B6u4V2rA|EY@_Dj(B
  integer, parameter :: im = 20         ! $B3J;RE@$N@_Dj(B (x)
  integer, parameter :: km = 20         ! $B3J;RE@$N@_Dj(B (z)

  ! $B:BI8JQ?t(B
  real(8), parameter :: xmin = 0.0D0
  real(8), parameter :: xmax = 1000.0D0
  real(8), parameter :: dx = (xmax-xmin)/im ! $B3J;R4V3V(B (x)
  real(8)            :: zs(0:im-1)       ! $BCO7A(B (z_s)
  real(8), parameter :: zmin = 0.0D0        ! $B%b%G%kNN0h$NDl(B
  real(8), parameter :: zmax = 1000.0D0     ! $B%b%G%kNN0h$N9b$5(B (H)
  real(8), parameter :: dz = (zmax-zmin)/km ! $B3J;R4V3V(B (z)

  !--- $B4pK\>l$NJQ?t(B
  real(8) :: bpi(0:im-1,0:km-1)            ! $BL5<!8505NO(B
  real(8) :: bpr(0:im-1,0:km-1)            ! $B05NO(B
  real(8) :: brh(0:im-1,0:km-1)            ! $BL)EY(B
  real(8) :: btm(0:im-1,0:km-1)         ! $B29EY(B
  real(8) :: bth(0:im-1,0:km-1)            ! $B290L(B
  real(8) :: bvt(0:im-1,0:km-1)            ! $B2>290L(B
  real(8) :: bqv(0:im-1,0:km-1)            ! $B?e>x5$$N:.9gHf(B

  !--- $B$=$NB>$NJQ?t(B
  real(8) :: g12(0:im-1,0:km-1)           ! $B%d%3%S%"%s$N(B1/2$B>h(B

  real(8) :: x(0:im-1,0:km-1)                        ! x $B:BI8(B
  real(8) :: z(0:im-1,0:km-1)                        ! z^{*} $B:BI8(B
  real(8) :: zn(0:im-1,0:km-1)                       ! z $B:BI8(B

  !--- $BJ*M}%Q%i%a!<%?(B
  real(8), parameter :: bpr0 = 1000.0D0    ! $B4pK\>l$NCOI=LL$G$N05NO(B
  real(8), parameter :: rd = 287.0D0       ! $B4%Ag6u5$$N5$BNDj?t(B
  real(8), parameter :: cp = 1952.0D0      ! $BDj05HfG.(B
  real(8), parameter :: cv = 1463.0D0      ! $BDj@QHfG.(B
  real(8), parameter :: grv = 9.8D0        ! $B=ENO(B

  !--- $B$=$NB>$NJQ?t(B
  integer            :: i, k                   ! do $BJQ?t(B

  !-------------------------------------------------
  ! X-Z $B:BI8$N@_Dj(B
  !-------------------------------------------------

  do k=0,km-1
     do i=0,im-1
        x(i,k) = dx*i
        !--- $BCO7A(B *
        zs(i) = 0.0D0
        zn(i,k) = dz*k
        z(i,k) = zmax*(zn(i,k) - zs(i))/(zmax - zs(i))
     end do
  end do

  call output_gtool4_init
  call output_gtool4

  !-------------------------------------------------
  ! $B4pK\>l$N7W;;(B
  !-------------------------------------------------

  !--- $BCOI=LL$NL5<!8505NO(B
  do i = 0,im-1
     bpr(i,0) = bpr0
!     bpr(i,0) = bpr0 - 0.005D0*x(i,0)
     bpi(i,0) = (bpr(i,0)/bpr0)**(rd/cp)
  end do

  !--- * $B$OCM$rM?$($k(B

  do k = 0,km-1
     do i = 0,im-1
        !--- $B29EY(B *
        btm(i,k) = 300.0D0 - 0.005D0*z(i,k)

        !--- $B?e>x5$$N:.9gHf(B *
        bqv(i,k) = 0.0D0

        !----------------------------------------------

        !--- $B%d%3%S%"%s$N(B1/2$B>h(B
        g12(i,k) = (zmax - zs(i))/zmax

        !--- $B05NO(B
        bpr(i,k) = bpr0*bpi(i,k)**(cp/rd)

        !--- $B290L(B
        bth(i,k) = btm(i,k)/bpi(i,k)
  
        !--- $B2>290L(B
        bvt(i,k) = bth(i,k)*(1.0D0 + 0.61D0*bqv(i,k))

        !--- $BL)EY(B
        brh(i,k) = (bpr0/rd)*(bpi(i,k)**(cv/cp))/bvt(i,k)

        call output_gtool4

     end do

     do i=0,im-1
        !--- $B@E?e05J?9U$N<0(B -> z+dz $B$G$NL5<!8505NO(B
        bpi(i,k+1) = bpi(i,k) - (z(i,k+1) - z(i,k))         &
                   & *g12(i,k)*grv/(cp*bvt(i,k))
     end do

  end do

  call output_gtool4_close

  stop

contains

!-----------------------------------------------------------
!--- $B%5%V%k!<%A%s(B

  !--- gtool4: $B%R%9%H%j!<=i4|2=(B
  subroutine output_gtool4_init
    call HistoryCreate(                                     & ! $B%R%9%H%j!<:n@.(B
          & file='base.nc',                                 &
          & title='2D cumulus model  base',                 &
          & source='kaminari - 2D cumulus model base',      &
          & institution='GFD_Dennou Club deepconv project', &
          & dims=(/'x','z'/), dimsizes=(/im,km/),           &
          & longnames=(/'X-coordinate', 'Z-coordinate'/),   &
          & units=(/'m','m'/), origin=0.0, interval=real(1.0)&
          & )

    call HistoryPut('x',x(0:im-1,0))                      ! $BJQ?t=PNO(B
    call HistoryPut('z',z(0,0:km-1))                      ! $BJQ?t=PNO(B

    call HistoryAddVariable(                            & ! $BJQ?tDj5A(B
         & varname='bpr', dims=(/'x','z'/),             & ! $B05NO(B
         & longname='pressure', units='hPa', xtype='double')
    call HistoryAddVariable(                            &
         & varname='btm', dims=(/'x','z'/),             & ! $B29EY(B
         & longname='temperature', units='K', xtype='double')
    call HistoryAddVariable(                            &
         & varname='brh', dims=(/'x','z'/),             & ! $BL)EY(B
         & longname='density', units='kg/m3',           &
         & xtype='double')
    call HistoryAddVariable(                            &
         & varname='bqv', dims=(/'x','z'/),             & ! $B?e>x5$:.9gHf(B
         & longname='mixing ratio of water vapor',      &
         &  units='g/g', xtype='double')
    call HistoryAddVariable(                            & ! $BJQ?tDj5A(B
         & varname='zs', dims=(/'x'/),                  & ! $BCO7A(B
         & longname='geometry', units='m', xtype='double')
  end subroutine output_gtool4_init

  !--- gtool4: $B%U%!%$%k=PNO(B
  subroutine output_gtool4
    write(6,*) 'i,k = ', i, k
    call HistoryPut('btm', btm)
    call HistoryPut('bpr', bpr)
    call HistoryPut('brh', brh)
    call HistoryPut('bqv', bqv)
    call HistoryPut('zs', zs)
  end subroutine output_gtool4

  !--- gtool4: $B=*N;=hM}(B
  subroutine output_gtool4_close
    call HistoryClose
  end subroutine output_gtool4_close

end program base
