c23456789012345678901234567890123456789012345678901234567890123456789012


      SUBROUTINE FFxxxxx 

c     ******************************************************************
c
c     This routine constructs the Jacobian matrix of the nonlinear 
c     functions
c
c     ******************************************************************

c     -------------
c     Declarations
c     -------------

c     ''''''''''''''''''''''
      implicit none
      include 'commonvar.f'      
c     ''''''''''''''''''''''

      integer	i,j,k,kcomp,n

      character*12 stringx

      character*12 stringr_swl,stringr_sel,stringr_seu,stringr_swu,
     >		   stringr_nwl,stringr_nel,stringr_neu,stringr_nwu
     
      character*12 stringU_w,stringU_e,
     >		   stringV_s,stringV_n,
     >		   stringW_l,stringW_u

      real*8	S,vol,fac




c     ----------------
c     Density equation
c     ----------------

c---- Initialization 

      kcomp = 0


c---- General loop on statements or grid boxes

c     ************
      do i=1,itot
      do j=1,jtot
      do k=1,ktot
c     ************ 


c---- Flag test

c++
      if (flag_den(i,j,k)) then
c++


c---- Statement counter

      kcomp = kcomp + 1


c---- Variables that appear in kcomp-th statement

c     ---------------------------------------
      stringr_swl = stringx('rr',i  ,j  ,k  )
      stringr_sel = stringx('rr',i+1,j  ,k  )
      stringr_seu = stringx('rr',i+1,j  ,k+1)
      stringr_swu = stringx('rr',i  ,j  ,k+1)
      stringr_nwl = stringx('rr',i  ,j+1,k  )
      stringr_nel = stringx('rr',i+1,j+1,k  )
      stringr_neu = stringx('rr',i+1,j+1,k+1)
      stringr_nwu = stringx('rr',i  ,j+1,k+1)
c--
      stringU_w = stringx('UU',i  ,j  ,k  )
      stringU_e = stringx('UU',i+1,j  ,k  )
      stringV_s = stringx('VV',i  ,j  ,k  )
      stringV_n = stringx('VV',i  ,j+1,k  )
      stringW_l = stringx('WW',i  ,j  ,k  )
      stringW_u = stringx('WW',i  ,j  ,k+1)
c     ---------------------------------------


c---- Elements of the kcomp-th row of FF

c++
      do n=1,n_o
c++

         if (string(n).eq.stringr_swl) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0
     >			  -U(i  ,j  ,k  ) / 4.d0
     >			  -W(i  ,j  ,k  ) / 4.d0

         else if (string(n).eq.stringr_sel) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0
     >                    +U(i+1,j  ,k  ) / 4.d0
     >                    -W(i  ,j  ,k  ) / 4.d0
 
         else if (string(n).eq.stringr_seu) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0
     >                    +U(i+1,j  ,k  ) / 4.d0
     >                    +W(i  ,j  ,k+1) / 4.d0

         else if (string(n).eq.stringr_swu) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0 
     >                    -U(i  ,j  ,k  ) / 4.d0
     >                    +W(i  ,j  ,k+1) / 4.d0 

         else if (string(n).eq.stringr_nwl) then

            FF(kcomp,n) = +V(i  ,j+1,k  ) / 4.d0 
     >                    -U(i  ,j  ,k  ) / 4.d0 
     >                    -W(i  ,j  ,k  ) / 4.d0 

         else if (string(n).eq.stringr_nel) then

            FF(kcomp,n) = +U(i+1,j  ,k  ) / 4.d0 
     >                    +V(i  ,j+1,k  ) / 4.d0 
     >                    -W(i  ,j  ,k  ) / 4.d0 

         else if (string(n).eq.stringr_neu) then

            FF(kcomp,n) = +U(i+1,j  ,k  ) / 4.d0 
     >                    +V(i  ,j+1,k  ) / 4.d0 
     >                    +W(i  ,j  ,k+1) / 4.d0 

         else if (string(n).eq.stringr_nwu) then

            FF(kcomp,n) = +V(i  ,j+1,k  ) / 4.d0 
     >                    -U(i  ,j  ,k  ) / 4.d0 
     >                    +W(i  ,j  ,k+1) / 4.d0 

         else if (string(n).eq.stringU_w) then

            FF(kcomp,n) = -(r(i  ,j  ,k  )+
     >                      r(i  ,j+1,k  )+
     >                      r(i  ,j+1,k+1)+
     >                      r(i  ,j  ,k+1)) / 4.d0 

         else if (string(n).eq.stringU_e) then

            FF(kcomp,n) = +(r(i+1,j  ,k  )+
     >                      r(i+1,j+1,k  )+
     >                      r(i+1,j+1,k+1)+
     >                      r(i+1,j  ,k+1)) / 4.d0 

         else if (string(n).eq.stringV_s) then

            FF(kcomp,n) = -(r(i  ,j  ,k  )+
     >                      r(i+1,j  ,k  )+
     >                      r(i+1,j  ,k+1)+
     >                      r(i  ,j  ,k+1)) / 4.d0 

         else if (string(n).eq.stringV_n) then

            FF(kcomp,n) = +(r(i  ,j+1,k  )+
     >                      r(i+1,j+1,k  )+
     >                      r(i+1,j+1,k+1)+
     >                      r(i  ,j+1,k+1)) / 4.d0 

         else if (string(n).eq.stringW_l) then

            FF(kcomp,n) = -(r(i  ,j  ,k  )+
     >                      r(i+1,j  ,k  )+
     >                      r(i+1,j+1,k  )+
     >                      r(i  ,j+1,k  )) / 4.d0 

         else if (string(n).eq.stringW_u) then

            FF(kcomp,n) = +(r(i  ,j  ,k+1)+
     >                      r(i+1,j  ,k+1)+
     >                      r(i+1,j+1,k+1)+
     >                      r(i  ,j+1,k+1)) / 4.d0 

         else 

            FF(kcomp,n) = 0.d0

         endif

c++
      enddo
c++


c---- End of flag test

c++
      endif
c++


c---- End of general loop

c     ***** 
      enddo
      enddo
      enddo
c     ***** 




c     ---------------
c     Tracer equation
c     ---------------

c---- General loop on statements or grid boxes

c     ************
      do i=1,itot
      do j=1,jtot
      do k=1,ktot
c     ************ 


c---- Flag test

c++
      if (flag_pos(i,j,k)) then
c++


c---- Factor

      S   = (r_0**2) * (sin(lat(j+1))-sin(lat(j))) * dx(i)
      vol = S * dz(k) 
      fac = prod * vol * sv 


c---- Statement counter

      kcomp = kcomp + 1


c---- Variables that appear in kcomp-th statement

c     ---------------------------------------
      stringr_swl = stringx('pp',i  ,j  ,k  )
      stringr_sel = stringx('pp',i+1,j  ,k  )
      stringr_seu = stringx('pp',i+1,j  ,k+1)
      stringr_swu = stringx('pp',i  ,j  ,k+1)
      stringr_nwl = stringx('pp',i  ,j+1,k  )
      stringr_nel = stringx('pp',i+1,j+1,k  )
      stringr_neu = stringx('pp',i+1,j+1,k+1)
      stringr_nwu = stringx('pp',i  ,j+1,k+1)
c--
      stringU_w = stringx('UU',i  ,j  ,k  )
      stringU_e = stringx('UU',i+1,j  ,k  )
      stringV_s = stringx('VV',i  ,j  ,k  )
      stringV_n = stringx('VV',i  ,j+1,k  )
      stringW_l = stringx('WW',i  ,j  ,k  )
      stringW_u = stringx('WW',i  ,j  ,k+1)
c     ---------------------------------------


c---- Elements of the kcomp-th row of FF

c++
      do n=1,n_o
c++

         if (string(n).eq.stringr_swl) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0
     >                    -U(i  ,j  ,k  ) / 4.d0
     >                    -W(i  ,j  ,k  ) / 4.d0  +  fac / 8.d0 

         else if (string(n).eq.stringr_sel) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0
     >                    +U(i+1,j  ,k  ) / 4.d0
     >                    -W(i  ,j  ,k  ) / 4.d0  +  fac / 8.d0
 
         else if (string(n).eq.stringr_seu) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0
     >                    +U(i+1,j  ,k  ) / 4.d0
     >                    +W(i  ,j  ,k+1) / 4.d0  +  fac / 8.d0

         else if (string(n).eq.stringr_swu) then

            FF(kcomp,n) = -V(i  ,j  ,k  ) / 4.d0 
     >                    -U(i  ,j  ,k  ) / 4.d0
     >                    +W(i  ,j  ,k+1) / 4.d0  +  fac / 8.d0

         else if (string(n).eq.stringr_nwl) then

            FF(kcomp,n) = +V(i  ,j+1,k  ) / 4.d0 
     >                    -U(i  ,j  ,k  ) / 4.d0 
     >                    -W(i  ,j  ,k  ) / 4.d0  +  fac / 8.d0 

         else if (string(n).eq.stringr_nel) then

            FF(kcomp,n) = +U(i+1,j  ,k  ) / 4.d0 
     >                    +V(i  ,j+1,k  ) / 4.d0 
     >                    -W(i  ,j  ,k  ) / 4.d0  +  fac / 8.d0 

         else if (string(n).eq.stringr_neu) then

            FF(kcomp,n) = +U(i+1,j  ,k  ) / 4.d0 
     >                    +V(i  ,j+1,k  ) / 4.d0 
     >                    +W(i  ,j  ,k+1) / 4.d0  +  fac / 8.d0 

         else if (string(n).eq.stringr_nwu) then

            FF(kcomp,n) = +V(i  ,j+1,k  ) / 4.d0 
     >                    -U(i  ,j  ,k  ) / 4.d0 
     >                    +W(i  ,j  ,k+1) / 4.d0  +  fac / 8.d0 

         else if (string(n).eq.stringU_w) then

            FF(kcomp,n) = -(p(i  ,j  ,k  )+
     >                      p(i  ,j+1,k  )+
     >                      p(i  ,j+1,k+1)+
     >                      p(i  ,j  ,k+1)) / 4.d0 

         else if (string(n).eq.stringU_e) then

            FF(kcomp,n) = +(p(i+1,j  ,k  )+
     >                      p(i+1,j+1,k  )+
     >                      p(i+1,j+1,k+1)+
     >                      p(i+1,j  ,k+1)) / 4.d0 

         else if (string(n).eq.stringV_s) then

            FF(kcomp,n) = -(p(i  ,j  ,k  )+
     >                      p(i+1,j  ,k  )+
     >                      p(i+1,j  ,k+1)+
     >                      p(i  ,j  ,k+1)) / 4.d0 

         else if (string(n).eq.stringV_n) then

            FF(kcomp,n) = +(p(i  ,j+1,k  )+
     >                      p(i+1,j+1,k  )+
     >                      p(i+1,j+1,k+1)+
     >                      p(i  ,j+1,k+1)) / 4.d0 

         else if (string(n).eq.stringW_l) then

            FF(kcomp,n) = -(p(i  ,j  ,k  )+
     >                      p(i+1,j  ,k  )+
     >                      p(i+1,j+1,k  )+
     >                      p(i  ,j+1,k  )) / 4.d0 

         else if (string(n).eq.stringW_u) then

            FF(kcomp,n) = +(p(i  ,j  ,k+1)+
     >                      p(i+1,j  ,k+1)+
     >                      p(i+1,j+1,k+1)+
     >                      p(i  ,j+1,k+1)) / 4.d0 

         else 

            FF(kcomp,n) = 0.d0

         endif

c++
      enddo
c++


c---- End of flag test

c++
      endif
c++


c---- End of general loop

c     ***** 
      enddo
      enddo
      enddo
c     ***** 




c     --------------------------------
c     Testing the number of rows of F 
c     --------------------------------

      if (kcomp.ne.n_nl) then
         write(*,*)
         write(*,*) 'f and Cf are not conformable'
         write(*,*)
         stop 'STOP'
      endif




c     ----------------------
c     End of the subroutine
c     ----------------------

      return

      end

c     ******************************************************************
