! This windset macro provides sky and north, south, east, west
   ! boundary conditions, based on nx, ny, nz and angwind

mesg(windset macro has been called   
  GROUP 7
  save7begin
solve(w1,p1)
solutn(p1,y,y,y,p,p,p)            
solutn(w1,y,y,y,p,p,p)            

if(nx.gt.1) then
solveu=t
solve(u1)
solutn(u1,y,y,y,p,p,p)  
endif          

if(ny.gt.1) then
solvev=t
solve(v1)
solutn(v1,y,y,y,p,p,p)  
endif                                             
  reference adiabatic-atmosphere settings
  ***************************************
  if(:atmos:.eq.adiabat) then
t0                            ! temperature
const1
const1=8314.0*t0/(29*9.81)    ! z0 = 'isothermal scale height'
const1
const2=1.4/(1.4-1.0)          ! 1.4 = gamma = specific-heat ratio
const1=1.0/(const1*const2)    ! multiplier of height = (gamma-1)/(gamma*z0)
 t0
 form=:t0:*(1-zg*:const1: )   ! formula fr reference temperaturs
(stored var tref is :form:)  
 p0
 form=:p0:*(1-zg*:const1:)^:const2: ! formula fr reference temperaturs
(stored var pref is :form:)
 d0=p0*29/(8314.0*t0)
 d0
 form=d0*(1-zg*:const1:)  
  (stored var dref is :form:)  
  endif  
  save7end

 angwind=angwind*pi/180                 ! u and v components
 uwind=velwind*cos(angwind)
 vwind=velwind*sin(angwind)
 uwind
 vwind

twonz=2*nz                  ! increase zwlast so the the velocity 
                            ! values to be easily checked
twonz
zwlast=zwlast*(twonz+1)/twonz
zwlast
           
  GROUP 11          
if(nx.gt.1) then
  fiinit(u1)=uwind
 (initial of u1 is uwind*zg/:zwlast:*:twonz:/(:twonz:-1))
endif  
if(ny.gt.1) then
  fiinit(v1)=vwind
 (initial of v1 is vwind*zg/:zwlast:*:twonz:/(:twonz:-1))
endif
coef=enut*twonz/zwlast ! This value corresponds to the enut set 
           ! above so as to provide a linear profile
           ! for testing the method                                 
coef2=coef*0.6667 ! But near boundaries one must multiply by 0.6667
                  ! because the diffusion 
                  ! term takes the area to be proportional to
                  ! dxu everywhere, whereas the source term has
                  ! uses an extrapolation to the domain boundary.
coef
coef2
   --------------------------------------- High (sky) boundary
patch(highp,high,1,nx,1,ny,nz,nz,1,lstep)  
coval(highp,p1,fixval,0.0)           ! this fixes the pressure      
                                     ! and horizontal          
                                     ! velocities in the       
                                     ! top (sky) slab          
if(solveu) then                          ! for u1
 coval(highp,u1,fixval,uwind) 
                
 patch(hiCEN_u,high,2,nx-2,1,ny,nz,nz,1,1)      ! central region
 coval(hiCEN_u,u1, fixval, uwind)            

 if(nx.gt.1) then
 patch(hiW_u,high,1,1,1,ny,nz,nz,1,1)        ! west edge
 coval(hiW_u,u1, fixval*factor, uwind)            
 
 patch(hiE_u,high,nx-1,nx,1,ny,nz,nz,1,1)  ! east edge
 coval(hiE_u,u1, fixval*factor, uwind)            
 endif
 
 if(ny.gt.1) then
 patch(hiN_u,high,2,nx-2,ny,ny,nz,nz,1,1)      ! north edge
 coval(hiN_u,u1,fixval*factor,uwind)
 
 patch(hiS_u,high,2,nx-2,1,1,nz,nz,1,1)        ! south edge
 coval(hiS_u,u1,fixval*factor,uwind)
 endif
endif

if(solvev) then                          ! for v1
 coval(highp,v1,fixval,vwind) 
 
 patch(hiCEN_v,high,1,nx,2,ny-2,nz,nz,1,1)      ! central region
 coval(hiCEN_v,v1, fixval, vwind)            

 if(ny.gt.1) then
 patch(hiS_v,high,1,nx,1,1,nz,nz,1,1)        ! south edge
 coval(hiS_v,v1, fixval*factor, vwind)            

 patch(hiN_v,high,1,nx,ny-1,ny,nz,nz,1,1)  ! north edge
 coval(hiN_v,v1, fixval*factor, vwind)            
 endif

 if(nx.gt.1) then
 patch(hiE_v,high,nx,nx,2,ny-2,nz,nz,1,1)      ! east edge
 coval(hiE_v,v1,fixval*factor,vwind)
 
 patch(hiW_v,high,1,1,2,ny-2,nz,nz,1,1)       ! west edge
 coval(hiW_v,v1,fixval*factor,vwind)
 endif
endif

 patch(highpm1,high,1,nx,1,ny,nz-1,nz-1,1,lstep)  ! this may 
 coval(highpm1,w1,1.e5,same)            ! improve convergence

            --------------------------------------- east boundary
if(nx.gt.1) then 
 patch(eastp,east,nx,nx,1,ny,1,nz,1,lstep)
 coval(eastp,p1,fixval,0.0)
 coval(eastp,w1,fixval,0.0)
 
 patch(&dfH_Eu,cell,nx-1,nx,1,ny,1,nz-1,1,lstep)  ! dominant high-diff.
 coval(&dfH_Eu,u1,factor,0.0)                       ! for u1
 
  patch(&sor_Eu,cell,nx-1,nx-1,1,ny,1,nz-1,1,lstep)  ! no source
  coval(&sor_Eu,u1,0.0,0.0)                          ! for u1
 
if(ny.gt.1) then 
 patch(&dfH_Ev,cell,nx,nx,1,ny,1,nz-1,1,lstep)  ! dominant high-diff.
 coval(&dfH_Ev,v1,factor,0.0)                   ! for v1
endif  
endif
                                ! coriolis source ( to be worked on )
                          
  coval(&sor_Eu,u1,factor,0.0)  ! note that this source-
                                ! multiplication had a bad effect on 
                                ! convergence because it increased 
                                ! the dveldps because of a line in 
                                ! chapter 9.6 of comp, the purpose 
                                ! of which is not yet clear.
  patch(&sor_Ev,cell,nx,nx,1,ny,1,nz-1,1,lstep)   ! multiply coriolis
  
  coval(&sor_Ev,v1,factor,0.0)
  if(itmod.eq.1) then
  coval(&dfH_E,ke,factor,0.0)
  coval(&dfH_E,ep,factor,0.0)
  endif


            --------------------------------------- North boundary
if(ny.gt.1) then 
 patch(northp,north,1,nx,ny,ny,1,nz,1,lstep)
 coval(northp,p1,fixval,0.0)
 
 patch(&dfH_Nv,cell,1,nx,ny-1,ny,1,nz-1,1,lstep)  ! dominant high-diff.
 coval(&dfH_Nv,v1,factor,0.0)                   ! for v1
if(nx.gt.1) then
 patch(&dfH_Nu,cell,1,nx,ny,ny,1,nz-1,1,lstep)  ! dominant high-diff.
 coval(&dfH_Nu,u1,factor,0.0)                       ! for u1
endif 
endif  
  
            --------------------------------------- 
            West boundary
if(nx.gt.1) then
 patch(westp,west,1,1,1,ny,1,nz,1,lstep)      ! pressure           
 coval(westp,p1,fixval,0.0)                                        
 coval(westp,w1,fixval,0.0)                                        
                                                                  
 patch(&dfH_W,cell,1,1,1,ny,1,nz-1,1,lstep)   ! dominant high-diff.
 coval(&dfH_W,u1,factor,0.0)                  ! for u1
 
  patch(&sor_W,cell,1,1,1,ny,1,nz-1,1,lstep)   ! no source for u1
  coval(&sor_W,u1,0.0,0.0)                  
 
 if(solvev) then
 coval(&dfH_W,v1,factor,0.0)        ! dominant high-diff. for v1
  patch(&sor_v1,cell,1,nx,1,1,1,nz-1,1,lstep)   ! no source
  coval(&sor_v1,v1,0.0,0.0)                     ! for v1
 endif
endif

            --------------------------------------- 
            South boundary
if(ny.gt.1) then
 patch(southp,south,1,nx,1,1,1,nz,1,lstep)    ! pressure           
 coval(southp,p1,fixval,0.0)                                        
 coval(southp,w1,fixval,0.0)                                        
                                                                  
 patch(&dfH_S,cell,1,nx,1,1,1,nz-1,1,lstep)   ! dominant high-diff.
 coval(&dfH_S,v1,factor,0.0)                  ! for v1
 
  patch(&sor_S,cell,1,nx,1,1,1,nz-1,1,lstep)   ! no source
  coval(&sor_S,v1,0.0,0.0)                     ! for v1
 
 if(solveu) then
 coval(&dfH_S,u1,factor,0.0)         ! dominant high diff. for u1
  coval(&sor_S,u1,0.0,0.0)            ! no source for u1
 endif
endif

            --------------------------------------- Low boundary
if(solveu) then                          ! for u1
 patch(lowCEN_u,low,2,nx-2,1,ny,1,1,1,1)      ! central region
 coval(lowCEN_u,u1, coef, 0.0)            

 if(nx.gt.1) then
 patch(lowW_u,low,1,1,1,ny,1,1,1,1)        ! west edge
 coval(lowW_u,u1, coef2*factor, 0.0)            
 
 patch(lowE_u,low,nx-1,nx,1,ny,1,1,1,1)  ! east edge
 coval(lowE_u,u1, coef2*factor, 0.0)            
 endif
 
 if(ny.gt.1) then
 patch(lowN_u,low,2,nx-2,ny,ny,1,1,1,1)      ! north edge
 coval(lowN_u,u1,coef*factor,0.0)
 
 patch(lowS_u,low,2,nx-2,1,1,1,1,1,1)        ! south edge
 coval(lowS_u,u1,coef*factor,0.0)
 endif
endif

if(solvev) then                          ! for v1
 patch(lowCEN_v,low,1,nx,2,ny-2,1,1,1,1)      ! central region
 coval(lowCEN_v,v1, coef, 0.0)            

 if(ny.gt.1) then
 patch(lowS_v,low,1,nx,1,1,1,1,1,1)        ! south edge
 coval(lowS_v,v1, coef2*factor, 0.0)            

 patch(lowN_v,low,1,nx,ny-1,ny,1,1,1,1)  ! north edge
 coval(lowN_v,v1, coef2*factor, 0.0)            
 endif

 if(nx.gt.1) then
 patch(lowE_v,low,nx,nx,2,ny-2,1,1,1,1)      ! east edge
 coval(lowE_v,v1,coef*factor,0.0)
 
 patch(lowW_v,low,1,1,2,ny-2,1,1,1,1)       ! west edge
 coval(lowW_v,v1,coef*factor,0.0)
 endif
endif