copyright (C) 2001 MSC-RPN COMM %%%MC2%%%
***s/r mtn_3d
*

      subroutine mtn_3d 1,9
      implicit none
*
#include "nesting.cdk"
#include "dynmem.cdk"
#include "tracers.cdk"
#include "nestpnt.cdk"
#include "yomdyn.cdk"
#include "yomdyn1.cdk"
#include "dynpar.cdk"
#include "consdyn_8.cdk"
#include "levels.cdk"
#include "partopo.cdk"
#include "mtn.cdk"
#include "topo.cdk"
#include "refer.cdk"
#include "h_geop.cdk"
#include "bcsdim.cdk"
#include "bcsmem.cdk"
*    
      integer i,j,k,n,dim,kc,kp1,km1,
     $        id,indn,indw,inde,ofn,ofw,ofe
      real s,d
      pointer (pas1 ,s(*)), (pad1 ,d(*))
      real*8 p00,one,dh,dqdz,pt5,a00,a01,ttt
*
*---------------------------------------------------------------------
*
      call hauteurs
*
      one  = 1.0d0
      pt5  = 0.5d0
*
* define the initial conditions for the winds
*
      uup = mtn_flo
      vvp = 0.
*
* define the initial conditions for temperature and pressure
*
      if (theocase.eq.'MTN_PINTY1'.or.theocase.eq.'MTN_PINTY2') then
*
*  pinty cases
* 
         bbp = 0.
         ppp = 0.
*
      else if (theocase.eq.'MTN_SHAER') then
*
*  shaer case
*
         tzero1 = 288.15
         nstar1 = 1.e-2
         p00    = 100000.0d0
         a00    = nstar1**2/grav_8
         a01    = grav_8/(cpd_8*tzero1)/a00
*
         call qntstar(qstr,nssq,gots,orts,ht,hm,
     $                 (maxx-minx+1)*(maxy-miny+1),0,gnk)
*
         do k=0,gnk
            do j=1-hy,ldnj+hy
            do i=1-hx,ldni+hx
               ppp(i,j,k) = dlog(p00) + cpd_8/rgasd_8*
     $                   dlog((one-a01)+a01*exp(-a00*hm(i,j,k)))
               ppp(i,j,k) = rgasd_8*grtstar*
     $                               (ppp(i,j,k) - qstr(1,1,k))
            end do
            end do
         end do
*
         do k=1,gnk
            km1=max(k-1,1)
            kp1=min(k+1,gnk)
            do j=1-hy,ldnj+hy
            do i=1-hx,ldni+hx
               dh=(hw(i,j,kp1)-hw(i,j,km1))*pt5
               dqdz=(ppp(i,j,k)-ppp(i,j,k-1))/dh
               bbp(i,j,k) = grav_8*dqdz/(grav_8-dqdz)
            end do
            end do
         end do
*
      else if (theocase.eq.'NOFLOW') then
*
* is it noflow?
*
         tzero1 = 175.
         tzero1 = 375.
         p00    = 100000.0d0
*
         do j=1-hy,ldnj+hy
         do i=1-hx,ldni+hx
            ppp(i,j,0) = 0.0
         end do
         end do
*
         do k=1,gnk
            km1=max(k-1,1)
            kp1=min(k+1,gnk)
            do j=1-hy,ldnj+hy
            do i=1-hx,ldni+hx
               bbp(i,j,k) = grav_8*(tzero1-grtstar)/grtstar
               ppp(i,j,k) = hm(i,j,k)*grav_8*(tzero1-grtstar)/tzero1
            end do
            end do
         end do
*
      endif
*
* Initialize Moisture and Vertical motion to zero and tracers to 1
*
      hup = 0.
      wwp = 0.
      trp = 1.
*
      dim = ndynvar*dim3d+dim2d
      pas1 = pappp
      pad1 = papp0
      do i=1,dim
         d(i) = s(i)
      end do
*
      ofn  = (bcs_in-bcs_is)/gnk
      ofw  = (bcs_iw-bcs_in)/gnk
      ofe  = (bcs_ie-bcs_iw)/gnk
      indn = bcs_in + ofn
      indw = bcs_iw + ofn + ofw
      inde = bcs_ie + ofn + ofw + ofe
      call trnes (uup,bcs_uu(bcs_is),bcs_uu(bcs_in),bcs_uu(bcs_iw),
     $                          bcs_uu(bcs_ie),minx,maxx,miny,maxy,
     $         minxs,maxxs,minys,maxys,minxw,maxxw,minyw,maxyw,gnk)
      call trnes (vvp,bcs_vv(bcs_is),bcs_vv(bcs_in),bcs_vv(bcs_iw),
     $                          bcs_vv(bcs_ie),minx,maxx,miny,maxy,
     $         minxs,maxxs,minys,maxys,minxw,maxxw,minyw,maxyw,gnk)
      call trnes (wwp,bcs_ww(bcs_is),bcs_ww(bcs_in),bcs_ww(bcs_iw),
     $                          bcs_ww(bcs_ie),minx,maxx,miny,maxy,
     $         minxs,maxxs,minys,maxys,minxw,maxxw,minyw,maxyw,gnk)
      call trnes (bbp,bcs_bb(bcs_is),bcs_bb(bcs_in),bcs_bb(bcs_iw),
     $                          bcs_bb(bcs_ie),minx,maxx,miny,maxy,
     $         minxs,maxxs,minys,maxys,minxw,maxxw,minyw,maxyw,gnk)
      call trnes (ppp,bcs_pp(bcs_is),bcs_pp(indn)  ,bcs_pp(indw)  ,
     $                          bcs_pp(inde)  ,minx,maxx,miny,maxy,
     $       minxs,maxxs,minys,maxys,minxw,maxxw,minyw,maxyw,gnk+1)
      call trnes (hup,bcs_hu(bcs_is),bcs_hu(bcs_in),bcs_hu(bcs_iw),
     $                          bcs_hu(bcs_ie),minx,maxx,miny,maxy,
     $         minxs,maxxs,minys,maxys,minxw,maxxw,minyw,maxyw,gnk)
      do n = 1, ntr
         id = (n-1)*bcs_sz+1
         call trnes (trp(1-hx,1-hy,1,n),bcs_tr(id),bcs_tr(id+bcs_in-1),
     $     bcs_tr(id+bcs_iw-1),bcs_tr(id+bcs_ie-1),minx,maxx,miny,maxy,
     $             minxs,maxxs,minys,maxys,minxw,maxxw,minyw,maxyw,gnk)
      end do
      if (ctebcs) then
         do i=1,(ntr+6)*bcs_sz+bcs_sz/gnk
            bcs_ppa(i) = bcs_pp(i)
         end do
      endif
*
      if (gnpilver.gt.0) then
      do k=1,gnk
      do j=1-hy,ldnj+hy
      do i=1-hx,ldni+hx
         ppntt(i,j,k) = ppp(i,j,k)
         uuntt(i,j,k) = uup(i,j,k)
         vvntt(i,j,k) = vvp(i,j,k)
         wwntt(i,j,k) = wwp(i,j,k)
         bbntt(i,j,k) = bbp(i,j,k)
         huntt(i,j,k) = hup(i,j,k)
      end do
      end do
      end do
      do n = 1, ntr
         do k=1,gnk
         do j=1-hy,ldnj+hy
         do i=1-hx,ldni+hx
            trntt(i,j,k,n) = trp(i,j,k,n)
         end do
         end do
         end do
      end do
      do j=1-hy,ldnj+hy
      do i=1-hx,ldni+hx
         ppntt(i,j,0) = ppp(i,j,0)
      end do
      end do
*
      if (ctebcs) then
         do k=1,gnk
         do j=1-hy,ldnj+hy
         do i=1-hx,ldni+hx
            ppnta(i,j,k) = ppp(i,j,k)
            uunta(i,j,k) = uup(i,j,k)
            vvnta(i,j,k) = vvp(i,j,k)
            wwnta(i,j,k) = wwp(i,j,k)
            bbnta(i,j,k) = bbp(i,j,k)
            hunta(i,j,k) = hup(i,j,k)
         end do
         end do
         do n=1,ntr
            do j=1-hy,ldnj+hy
            do i=1-hx,ldni+hx
               trnta(i,j,k,n) = trp(i,j,k,n)
            end do
            end do
         end do
         end do
         do j=1-hy,ldnj+hy
         do i=1-hx,ldni+hx
            ppnta(i,j,0) = ppp(i,j,0)
         end do
         end do
      endif
      endif
*     
      if (myproc.eq.0) then
         write (6,602)
      endif
*
 602  format (/'INITIALIZATION OF MOUNTAIN PROBLEM COMPLETED')
*---------------------------------------------------------------------
      return
      end