! ======================================================================================
! This file was generated by the version 5.3.6 of DFT on 08/15/2010. The differentiation
! transforming system(DFT) was jointly developed and sponsored by LASG of IAP(1998-2010)
! and LSEC of ICMSEC, AMSS(2001-2003)
! The copyright of the DFT system was declared by Walls at LASG, 1998-2010
! ======================================================================================

 MODULE g_module_sfs_nba

 USE module_configure, ONLY : grid_config_rec_type
 USE module_sfs_nba, ONLY : c1, c2, c3, ce, cb, cs  ! Added by Ning Pan, 2010-08-18

 IMPLICIT NONE

! REAL :: c1,c2,c3,ce,cb,cs  ! Remarked by Ning Pan, 2010-08-18

 CONTAINS

! Remarked by Ning Pan, 2010-08-18
! SUBROUTINE g_calc_mij_constants()

! IMPLICIT NONE

! REAL :: Tmpv1,g_Tmpv1
! REAL :: sk,pi

! sk =0.5

! pi =3.1415927

! cb =0.36

! cs =((8.0 *(1.0+cb))/(27.0 *pi**2))**0.5

! c1 =((960.0**0.5) *cb)/(7.0 *(1.0+cb) *sk)

! c2 =c1

! ce =((8.0 *pi/27.0)**(1.0/3.0)) *cs**(4.0/3.0)

! c3 =((27.0/(8.0 *pi))**(1.0/3.0)) *cs**(2.0/3.0)
! Return

! END SUBROUTINE g_calc_mij_constants

 SUBROUTINE g_calc_smnsmn(smnsmn,g_smnsmn,s11,g_s11,s22,g_s22,s33, &
 g_s33,s12,g_s12,s13,g_s13,s23,g_s23,config_flags,ids,ide,jds,jde,kds,kde, &
 ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,its,ite,jts,jte,kts,kte)

 IMPLICIT NONE

 REAL :: Tmpv1,g_Tmpv1,Tmpv2,g_Tmpv2,Tmpv3,g_Tmpv3
 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: smnsmn,g_smnsmn

 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: s11,g_s11,s22,g_s22,s33,g_s33, &
 s12,g_s12,s13,g_s13,s23,g_s23
 TYPE(grid_config_rec_type) :: config_flags
 INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe, &
 its,ite,jts,jte,kts,kte
 REAL :: tmp,g_tmp
 INTEGER :: i,j,k,i_start,i_end,j_start,j_end,ktf

 ktf =min(kte,kde-1)

 i_start =its

 i_end =min(ite,ide-1)

 j_start =jts

 j_end =min(jte,jde-1)

 IF( config_flags%open_xs .or. config_flags%specified .or.   &
       config_flags%nested) i_start =max(ids+1,its)

 IF( config_flags%open_xe .or. config_flags%specified .or.   &
       config_flags%nested) i_end =min(ide-2,ite)

 IF( config_flags%open_ys .or. config_flags%specified .or.   &
       config_flags%nested) j_start =max(jds+1,jts)

 IF( config_flags%open_ye .or. config_flags%specified .or.   &
       config_flags%nested) j_end =min(jde-2,jte)

 IF( config_flags%periodic_x ) i_start =its

 IF( config_flags%periodic_x ) i_end =min(ite,ide-1)

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_Tmpv1 =2.0*s11(i,k,j)*g_s11(i,k,j) 
 Tmpv1 =s11(i,k,j)*s11(i,k,j)

 g_Tmpv2 =2.0*s22(i,k,j)*g_s22(i,k,j) 
 Tmpv2 =s22(i,k,j)*s22(i,k,j)

 g_Tmpv3 =2.0*s33(i,k,j)*g_s33(i,k,j) 
 Tmpv3 =s33(i,k,j)*s33(i,k,j)

 g_smnsmn(i,k,j) =0.25*(g_Tmpv1 +g_Tmpv2 +g_Tmpv3)
 smnsmn(i,k,j) =0.25*(Tmpv1 +Tmpv2 +Tmpv3)

 ENDDO
 ENDDO
 ENDDO

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_tmp =0.125*(g_s12(i,k,j) +g_s12(i,k,j+1) +g_s12(i+1,k,j) +g_s12(i+1,k,j+1))
 tmp =0.125*(s12(i,k,j) +s12(i,k,j+1) +s12(i+1,k,j) +s12(i+1,k,j+1))

 g_Tmpv1 =2.0*tmp*g_tmp +2.0*g_tmp*tmp 
 Tmpv1 =2.0*tmp*tmp

 g_smnsmn(i,k,j) =g_smnsmn(i,k,j) +g_Tmpv1
 smnsmn(i,k,j) =smnsmn(i,k,j) +Tmpv1

 ENDDO
 ENDDO
 ENDDO

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_tmp =0.125*(g_s13(i,k+1,j) +g_s13(i,k,j) +g_s13(i+1,k+1,j) +g_s13(i+1,k,j))
 tmp =0.125*(s13(i,k+1,j) +s13(i,k,j) +s13(i+1,k+1,j) +s13(i+1,k,j))

 g_Tmpv1 =2.0*tmp*g_tmp +2.0*g_tmp*tmp 
 Tmpv1 =2.0*tmp*tmp

 g_smnsmn(i,k,j) =g_smnsmn(i,k,j) +g_Tmpv1
 smnsmn(i,k,j) =smnsmn(i,k,j) +Tmpv1

 ENDDO
 ENDDO
 ENDDO

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_tmp =0.125*(g_s23(i,k+1,j) +g_s23(i,k,j) +g_s23(i,k+1,j+1) +g_s23(i,k,j+1))
 tmp =0.125*(s23(i,k+1,j) +s23(i,k,j) +s23(i,k+1,j+1) +s23(i,k,j+1))

 g_Tmpv1 =2.0*tmp*g_tmp +2.0*g_tmp*tmp 
 Tmpv1 =2.0*tmp*tmp

 g_smnsmn(i,k,j) =g_smnsmn(i,k,j) +g_Tmpv1
 smnsmn(i,k,j) =smnsmn(i,k,j) +Tmpv1

 ENDDO
 ENDDO
 ENDDO

 Return

 END SUBROUTINE g_calc_smnsmn

 SUBROUTINE g_calc_mii(m11,g_m11,m22,g_m22,m33,g_m33,s11,g_s11,s22, &
 g_s22,s33,g_s33,s12,g_s12,s13,g_s13,s23,g_s23,r12,g_r12,r13, &
 g_r13,r23,g_r23,smnsmn,g_smnsmn,tke,g_tke,rdzw,g_rdzw,dx,dy, &
 config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe, &
 its,ite,jts,jte,kts,kte)

 IMPLICIT NONE

 REAL :: Tmpv1,g_Tmpv1,Tmpv2,g_Tmpv2,Tmpv3,g_Tmpv3,Tmpv4,g_Tmpv4,Tmpv5, &
 g_Tmpv5,Tmpv6,g_Tmpv6,Tmpv7,g_Tmpv7,Tmpv8,g_Tmpv8

   REAL :: g_Sqrt
 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: m11,g_m11,m22,g_m22,m33,g_m33

 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: s11,g_s11,s22,g_s22,s33,g_s33, &
 s12,g_s12,s13,g_s13,s23,g_s23,r12,g_r12,r13,g_r13,r23,g_r23,smnsmn, &
 g_smnsmn,tke,g_tke,rdzw,g_rdzw

 REAL :: dx,dy
 TYPE(grid_config_rec_type) :: config_flags
 INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe, &
 its,ite,jts,jte,kts,kte

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: ss11,g_ss11,ss22,g_ss22, &
 ss33,g_ss33,ss12,g_ss12,ss13,g_ss13,ss23,g_ss23,rr12,g_rr12,rr13, &
 g_rr13,rr23,g_rr23

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: ss12c,g_ss12c,rr12c,g_rr12c, &
 ss13c,g_ss13c,rr13c,g_rr13c,ss23c,g_ss23c,rr23c,g_rr23c
 REAL :: delta,g_delta,a,g_a,b,g_b
 INTEGER :: i,j,k,i_start,i_end,j_start,j_end,ktf,is_ext,js_ext



 ktf =min(kte,kde-1)

 i_start =its

 i_end =ite

 j_start =jts

 j_end =jte

 IF( config_flags%open_xs .OR. config_flags%specified .OR.   &
         config_flags%nested) i_start =max(ids+1,its)

 IF( config_flags%open_xe .OR. config_flags%specified .OR.   &
         config_flags%nested) i_end =min(ide-1,ite)

 IF( config_flags%open_ys .OR. config_flags%specified .OR.   &
         config_flags%nested) j_start =max(jds+1,jts)

 IF( config_flags%open_ye .OR. config_flags%specified .OR.   &
         config_flags%nested) j_end =min(jde-1,jte)

 IF( config_flags%periodic_x ) i_start =its

 IF( config_flags%periodic_x ) i_end =ite

 is_ext =1

 js_ext =1

 i_start =i_start-is_ext

 j_start =j_start-js_ext

 DO j =j_start,j_end+1
 DO k =kts,ktf
 DO i =i_start,i_end+1

 g_ss11(i,k,j) =g_s11(i,k,j)/2.0
 ss11(i,k,j) =s11(i,k,j)/2.0

 g_ss22(i,k,j) =g_s22(i,k,j)/2.0
 ss22(i,k,j) =s22(i,k,j)/2.0

 g_ss33(i,k,j) =g_s33(i,k,j)/2.0
 ss33(i,k,j) =s33(i,k,j)/2.0

 g_ss12(i,k,j) =g_s12(i,k,j)/2.0
 ss12(i,k,j) =s12(i,k,j)/2.0

 g_ss13(i,k,j) =g_s13(i,k,j)/2.0
 ss13(i,k,j) =s13(i,k,j)/2.0

 g_ss23(i,k,j) =g_s23(i,k,j)/2.0
 ss23(i,k,j) =s23(i,k,j)/2.0

 g_rr12(i,k,j) =g_r12(i,k,j)/2.0
 rr12(i,k,j) =r12(i,k,j)/2.0

 g_rr13(i,k,j) =g_r13(i,k,j)/2.0
 rr13(i,k,j) =r13(i,k,j)/2.0

 g_rr23(i,k,j) =g_r23(i,k,j)/2.0
 rr23(i,k,j) =r23(i,k,j)/2.0

 ENDDO
 ENDDO
 ENDDO

 DO j =j_start,j_end+1
 DO i =i_start,i_end+1

 g_ss13(i,kde,j) =0.0
 ss13(i,kde,j) =0.0

 g_ss23(i,kde,j) =0.0
 ss23(i,kde,j) =0.0

 g_rr13(i,kde,j) =0.0
 rr13(i,kde,j) =0.0

 g_rr23(i,kde,j) =0.0
 rr23(i,kde,j) =0.0

 ENDDO
 ENDDO

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_ss12c(i,k,j) =0.25*(g_ss12(i,k,j) +g_ss12(i,k,j+1) +g_ss12(i+1,k,j) &
 +g_ss12(i+1,k,j+1))
 ss12c(i,k,j) =0.25*(ss12(i,k,j) +ss12(i,k,j+1) +ss12(i+1,k,j) +ss12(i+1,k,j+1))

 g_rr12c(i,k,j) =0.25*(g_rr12(i,k,j) +g_rr12(i,k,j+1) +g_rr12(i+1,k,j) &
 +g_rr12(i+1,k,j+1))
 rr12c(i,k,j) =0.25*(rr12(i,k,j) +rr12(i,k,j+1) +rr12(i+1,k,j) +rr12(i+1,k,j+1))

 g_ss13c(i,k,j) =0.25*(g_ss13(i,k+1,j) +g_ss13(i,k,j) +g_ss13(i+1,k+1,j) &
 +g_ss13(i+1,k,j))
 ss13c(i,k,j) =0.25*(ss13(i,k+1,j) +ss13(i,k,j) +ss13(i+1,k+1,j) +ss13(i+1,k,j))

 g_rr13c(i,k,j) =0.25*(g_rr13(i,k+1,j) +g_rr13(i,k,j) +g_rr13(i+1,k+1,j) &
 +g_rr13(i+1,k,j))
 rr13c(i,k,j) =0.25*(rr13(i,k+1,j) +rr13(i,k,j) +rr13(i+1,k+1,j) +rr13(i+1,k,j))

 g_ss23c(i,k,j) =0.25*(g_ss23(i,k+1,j) +g_ss23(i,k,j) +g_ss23(i,k+1,j+1) &
 +g_ss23(i,k,j+1))
 ss23c(i,k,j) =0.25*(ss23(i,k+1,j) +ss23(i,k,j) +ss23(i,k+1,j+1) +ss23(i,k,j+1))

 g_rr23c(i,k,j) =0.25*(g_rr23(i,k+1,j) +g_rr23(i,k,j) +g_rr23(i,k+1,j+1) &
 +g_rr23(i,k,j+1))
 rr23c(i,k,j) =0.25*(rr23(i,k+1,j) +rr23(i,k,j) +rr23(i,k+1,j+1) +rr23(i,k,j+1))

 ENDDO
 ENDDO
 ENDDO

 IF( config_flags%sfs_opt .EQ. 1 ) THEN

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0*2.0*(cs*g_delta)*(cs*delta)
 a =-1.0*(cs*delta)**2

 g_Tmpv1 =2.0*sqrt(2.0*smnsmn(i,k,j))*g_ss11(i,k,j) +2.0*g_Sqrt(2.0* &
 g_smnsmn(i,k,j), 2.0*smnsmn(i,k,j))*ss11(i,k,j) 
 Tmpv1 =2.0*sqrt(2.0*smnsmn(i,k,j))*ss11(i,k,j)

 g_Tmpv2 =2.0*ss11(i,k,j)*g_ss11(i,k,j) 
 Tmpv2 =ss11(i,k,j)*ss11(i,k,j)

 g_Tmpv3 =2.0*ss12c(i,k,j)*g_ss12c(i,k,j) 
 Tmpv3 =ss12c(i,k,j)*ss12c(i,k,j)

 g_Tmpv4 =2.0*ss13c(i,k,j)*g_ss13c(i,k,j) 
 Tmpv4 =ss13c(i,k,j)*ss13c(i,k,j)

 g_Tmpv5 =ss12c(i,k,j)*g_rr12c(i,k,j) +g_ss12c(i,k,j)*rr12c(i,k,j) 
 Tmpv5 =ss12c(i,k,j)*rr12c(i,k,j)

 g_Tmpv6 =ss13c(i,k,j)*g_rr13c(i,k,j) +g_ss13c(i,k,j)*rr13c(i,k,j) 
 Tmpv6 =ss13c(i,k,j)*rr13c(i,k,j)

 g_Tmpv7 =a*(g_Tmpv1 +c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4 -(g_smnsmn(i,k, &
 j)/3.0)) +c2*(-2.0*(g_Tmpv5 +g_Tmpv6))) +g_a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 + &
 Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(-2.0*(Tmpv5 +Tmpv6))) 
 Tmpv7 =a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(-2.0*(Tmpv5 +Tmpv6)))

 g_m11(i,k,j) =g_Tmpv7
 m11(i,k,j) =Tmpv7

 g_Tmpv1 =2.0*sqrt(2.0*smnsmn(i,k,j))*g_ss22(i,k,j) +2.0*g_Sqrt(2.0* &
 g_smnsmn(i,k,j), 2.0*smnsmn(i,k,j))*ss22(i,k,j) 
 Tmpv1 =2.0*sqrt(2.0*smnsmn(i,k,j))*ss22(i,k,j)

 g_Tmpv2 =2.0*ss22(i,k,j)*g_ss22(i,k,j) 
 Tmpv2 =ss22(i,k,j)*ss22(i,k,j)

 g_Tmpv3 =2.0*ss12c(i,k,j)*g_ss12c(i,k,j) 
 Tmpv3 =ss12c(i,k,j)*ss12c(i,k,j)

 g_Tmpv4 =2.0*ss23c(i,k,j)*g_ss23c(i,k,j) 
 Tmpv4 =ss23c(i,k,j)*ss23c(i,k,j)

 g_Tmpv5 =ss12c(i,k,j)*g_rr12c(i,k,j) +g_ss12c(i,k,j)*rr12c(i,k,j) 
 Tmpv5 =ss12c(i,k,j)*rr12c(i,k,j)

 g_Tmpv6 =ss23c(i,k,j)*g_rr23c(i,k,j) +g_ss23c(i,k,j)*rr23c(i,k,j) 
 Tmpv6 =ss23c(i,k,j)*rr23c(i,k,j)

 g_Tmpv7 =a*(g_Tmpv1 +c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4 -(g_smnsmn(i,k, &
 j)/3.0)) +c2*(2.0*(g_Tmpv5 -g_Tmpv6))) +g_a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 + &
 Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(2.0*(Tmpv5 -Tmpv6))) 
 Tmpv7 =a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(2.0*(Tmpv5 -Tmpv6)))

 g_m22(i,k,j) =g_Tmpv7
 m22(i,k,j) =Tmpv7

 g_Tmpv1 =2.0*sqrt(2.0*smnsmn(i,k,j))*g_ss33(i,k,j) +2.0*g_Sqrt(2.0* &
 g_smnsmn(i,k,j), 2.0*smnsmn(i,k,j))*ss33(i,k,j) 
 Tmpv1 =2.0*sqrt(2.0*smnsmn(i,k,j))*ss33(i,k,j)

 g_Tmpv2 =2.0*ss33(i,k,j)*g_ss33(i,k,j) 
 Tmpv2 =ss33(i,k,j)*ss33(i,k,j)

 g_Tmpv3 =2.0*ss13c(i,k,j)*g_ss13c(i,k,j) 
 Tmpv3 =ss13c(i,k,j)*ss13c(i,k,j)

 g_Tmpv4 =2.0*ss23c(i,k,j)*g_ss23c(i,k,j) 
 Tmpv4 =ss23c(i,k,j)*ss23c(i,k,j)

 g_Tmpv5 =ss13c(i,k,j)*g_rr13c(i,k,j) +g_ss13c(i,k,j)*rr13c(i,k,j) 
 Tmpv5 =ss13c(i,k,j)*rr13c(i,k,j)

 g_Tmpv6 =ss23c(i,k,j)*g_rr23c(i,k,j) +g_ss23c(i,k,j)*rr23c(i,k,j) 
 Tmpv6 =ss23c(i,k,j)*rr23c(i,k,j)

 g_Tmpv7 =a*(g_Tmpv1 +c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4 -(g_smnsmn(i,k, &
 j)/3.0)) +c2*(2.0*(g_Tmpv5 +g_Tmpv6))) +g_a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 + &
 Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(2.0*(Tmpv5 +Tmpv6))) 
 Tmpv7 =a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(2.0*(Tmpv5 +Tmpv6)))

 g_m33(i,k,j) =g_Tmpv7
 m33(i,k,j) =Tmpv7

 ENDDO
 ENDDO
 ENDDO

 ELSE

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0 *ce*g_delta
 a =-1.0 *ce*delta

 g_b =c3*g_delta
 b =c3*delta

 g_Tmpv1 =2.0*sqrt(tke(i,k,j))*g_ss11(i,k,j) +2.0*g_Sqrt(g_tke(i,k,j) &
, tke(i,k,j))*ss11(i,k,j) 
 Tmpv1 =2.0*sqrt(tke(i,k,j))*ss11(i,k,j)

 g_Tmpv2 =2.0*ss11(i,k,j)*g_ss11(i,k,j) 
 Tmpv2 =ss11(i,k,j)*ss11(i,k,j)

 g_Tmpv3 =2.0*ss12c(i,k,j)*g_ss12c(i,k,j) 
 Tmpv3 =ss12c(i,k,j)*ss12c(i,k,j)

 g_Tmpv4 =2.0*ss13c(i,k,j)*g_ss13c(i,k,j) 
 Tmpv4 =ss13c(i,k,j)*ss13c(i,k,j)

 g_Tmpv5 =ss12c(i,k,j)*g_rr12c(i,k,j) +g_ss12c(i,k,j)*rr12c(i,k,j) 
 Tmpv5 =ss12c(i,k,j)*rr12c(i,k,j)

 g_Tmpv6 =ss13c(i,k,j)*g_rr13c(i,k,j) +g_ss13c(i,k,j)*rr13c(i,k,j) 
 Tmpv6 =ss13c(i,k,j)*rr13c(i,k,j)

 g_Tmpv7 =b*(c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4 -(g_smnsmn(i,k,j)/3.0)) &
 +c2*(-2.0*(g_Tmpv5 +g_Tmpv6))) +g_b*(c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j) &
/3.0) +c2*(-2.0*(Tmpv5 +Tmpv6))) 
 Tmpv7 =b*(c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(-2.0*(Tmpv5 +Tmpv6)))

 g_Tmpv8 =a*(g_Tmpv1 +g_Tmpv7) +g_a*(Tmpv1 +Tmpv7) 
 Tmpv8 =a*(Tmpv1 +Tmpv7)

 g_m11(i,k,j) =g_Tmpv8
 m11(i,k,j) =Tmpv8

 g_Tmpv1 =2.0*sqrt(tke(i,k,j))*g_ss22(i,k,j) +2.0*g_Sqrt(g_tke(i,k,j) &
, tke(i,k,j))*ss22(i,k,j) 
 Tmpv1 =2.0*sqrt(tke(i,k,j))*ss22(i,k,j)

 g_Tmpv2 =2.0*ss22(i,k,j)*g_ss22(i,k,j) 
 Tmpv2 =ss22(i,k,j)*ss22(i,k,j)

 g_Tmpv3 =2.0*ss12c(i,k,j)*g_ss12c(i,k,j) 
 Tmpv3 =ss12c(i,k,j)*ss12c(i,k,j)

 g_Tmpv4 =2.0*ss23c(i,k,j)*g_ss23c(i,k,j) 
 Tmpv4 =ss23c(i,k,j)*ss23c(i,k,j)

 g_Tmpv5 =ss12c(i,k,j)*g_rr12c(i,k,j) +g_ss12c(i,k,j)*rr12c(i,k,j) 
 Tmpv5 =ss12c(i,k,j)*rr12c(i,k,j)

 g_Tmpv6 =ss23c(i,k,j)*g_rr23c(i,k,j) +g_ss23c(i,k,j)*rr23c(i,k,j) 
 Tmpv6 =ss23c(i,k,j)*rr23c(i,k,j)

 g_Tmpv7 =b*(c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4 -(g_smnsmn(i,k,j)/3.0)) &
 +c2*(2.0*(g_Tmpv5 -g_Tmpv6))) +g_b*(c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j) &
/3.0) +c2*(2.0*(Tmpv5 -Tmpv6))) 
 Tmpv7 =b*(c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(2.0*(Tmpv5 -Tmpv6)))

 g_Tmpv8 =a*(g_Tmpv1 +g_Tmpv7) +g_a*(Tmpv1 +Tmpv7) 
 Tmpv8 =a*(Tmpv1 +Tmpv7)

 g_m22(i,k,j) =g_Tmpv8
 m22(i,k,j) =Tmpv8

 g_Tmpv1 =2.0*sqrt(tke(i,k,j))*g_ss33(i,k,j) +2.0*g_Sqrt(g_tke(i,k,j) &
, tke(i,k,j))*ss33(i,k,j) 
 Tmpv1 =2.0*sqrt(tke(i,k,j))*ss33(i,k,j)

 g_Tmpv2 =2.0*ss33(i,k,j)*g_ss33(i,k,j) 
 Tmpv2 =ss33(i,k,j)*ss33(i,k,j)

 g_Tmpv3 =2.0*ss13c(i,k,j)*g_ss13c(i,k,j) 
 Tmpv3 =ss13c(i,k,j)*ss13c(i,k,j)

 g_Tmpv4 =2.0*ss23c(i,k,j)*g_ss23c(i,k,j) 
 Tmpv4 =ss23c(i,k,j)*ss23c(i,k,j)

 g_Tmpv5 =ss13c(i,k,j)*g_rr13c(i,k,j) +g_ss13c(i,k,j)*rr13c(i,k,j) 
 Tmpv5 =ss13c(i,k,j)*rr13c(i,k,j)

 g_Tmpv6 =ss23c(i,k,j)*g_rr23c(i,k,j) +g_ss23c(i,k,j)*rr23c(i,k,j) 
 Tmpv6 =ss23c(i,k,j)*rr23c(i,k,j)

 g_Tmpv7 =b*(c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4 -(g_smnsmn(i,k,j)/3.0)) &
 +c2*(2.0*(g_Tmpv5 +g_Tmpv6))) +g_b*(c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j) &
/3.0) +c2*(2.0*(Tmpv5 +Tmpv6))) 
 Tmpv7 =b*(c1*(Tmpv2 +Tmpv3 +Tmpv4 -smnsmn(i,k,j)/3.0) +c2*(2.0*(Tmpv5 +Tmpv6)))

 g_Tmpv8 =a*(g_Tmpv1 +g_Tmpv7) +g_a*(Tmpv1 +Tmpv7) 
 Tmpv8 =a*(Tmpv1 +Tmpv7)

 g_m33(i,k,j) =g_Tmpv8
 m33(i,k,j) =Tmpv8

 ENDDO
 ENDDO
 ENDDO
 ENDIF
 Return

 END SUBROUTINE g_calc_mii

 SUBROUTINE g_calc_m12(m12,g_m12,s11,g_s11,s22,g_s22,s12,g_s12,s13, &
 g_s13,s23,g_s23,r12,g_r12,r13,g_r13,r23,g_r23,smnsmn,g_smnsmn,tke, &
 g_tke,rdzw,g_rdzw,dx,dy,config_flags,ids,ide,jds,jde,kds,kde,ims,ime,jms,jme, &
 kms,kme,ips,ipe,jps,jpe,kps,kpe,its,ite,jts,jte,kts,kte)

 IMPLICIT NONE

 REAL :: Tmpv1,g_Tmpv1,Tmpv2,g_Tmpv2,Tmpv3,g_Tmpv3,Tmpv4,g_Tmpv4,Tmpv5, &
 g_Tmpv5,Tmpv6,g_Tmpv6,Tmpv7,g_Tmpv7,Tmpv8,g_Tmpv8,Tmpv9,g_Tmpv9, &
 Tmpv10,g_Tmpv10

 REAL :: g_Sqrt
 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: m12,g_m12

 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: s11,g_s11,s22,g_s22,s12,g_s12, &
 s13,g_s13,s23,g_s23,r12,g_r12,r13,g_r13,r23,g_r23,smnsmn,g_smnsmn, &
 tke,g_tke,rdzw,g_rdzw

 REAL :: dx,dy
 TYPE(grid_config_rec_type) :: config_flags
 INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe, &
 its,ite,jts,jte,kts,kte

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: ss11,g_ss11,ss22,g_ss22, &
 ss12,g_ss12,ss13,g_ss13,ss23,g_ss23,rr12,g_rr12,rr13,g_rr13,rr23,g_rr23

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: tked,g_tked,ss11d,g_ss11d, &
 ss22d,g_ss22d,ss13d,g_ss13d,ss23d,g_ss23d,rr13d,g_rr13d,rr23d,g_rr23d, &
 smnsmnd,g_smnsmnd
 REAL :: delta,g_delta,a,g_a,b,g_b
 INTEGER :: i,j,k,i_start,i_end,j_start,j_end,ktf,je_ext,ie_ext

 ktf =min(kte,kde-1)

 i_start =its

 i_end =ite

 j_start =jts

 j_end =jte

 IF( config_flags%open_xs .OR. config_flags%specified .OR.   &
         config_flags%nested ) i_start =max(ids+1,its)

 IF( config_flags%open_xe .OR. config_flags%specified .OR.   &
         config_flags%nested ) i_end =min(ide-1,ite)

 IF( config_flags%open_ys .OR. config_flags%specified .OR.   &
         config_flags%nested ) j_start =max(jds+1,jts)

 IF( config_flags%open_ye .OR. config_flags%specified .OR.   &
         config_flags%nested ) j_end =min(jde-1,jte)

 IF( config_flags%periodic_x ) i_start =its

 IF( config_flags%periodic_x ) i_end =ite

 je_ext =1

 ie_ext =1

 i_end =i_end+ie_ext

 j_end =j_end+je_ext

 DO j =j_start-1,j_end
 DO k =kts,ktf
 DO i =i_start-1,i_end

 g_ss11(i,k,j) =g_s11(i,k,j)/2.0
 ss11(i,k,j) =s11(i,k,j)/2.0

 g_ss22(i,k,j) =g_s22(i,k,j)/2.0
 ss22(i,k,j) =s22(i,k,j)/2.0

 g_ss12(i,k,j) =g_s12(i,k,j)/2.0
 ss12(i,k,j) =s12(i,k,j)/2.0

 g_ss13(i,k,j) =g_s13(i,k,j)/2.0
 ss13(i,k,j) =s13(i,k,j)/2.0

 g_ss23(i,k,j) =g_s23(i,k,j)/2.0
 ss23(i,k,j) =s23(i,k,j)/2.0

 g_rr12(i,k,j) =g_r12(i,k,j)/2.0
 rr12(i,k,j) =r12(i,k,j)/2.0

 g_rr13(i,k,j) =g_r13(i,k,j)/2.0
 rr13(i,k,j) =r13(i,k,j)/2.0

 g_rr23(i,k,j) =g_r23(i,k,j)/2.0
 rr23(i,k,j) =r23(i,k,j)/2.0

 ENDDO
 ENDDO
 ENDDO

 DO j =j_start-1,j_end
 DO i =i_start-1,i_end

 g_ss13(i,kde,j) =0.0
 ss13(i,kde,j) =0.0

 g_ss23(i,kde,j) =0.0
 ss23(i,kde,j) =0.0

 g_rr13(i,kde,j) =0.0
 rr13(i,kde,j) =0.0

 g_rr23(i,kde,j) =0.0
 rr23(i,kde,j) =0.0

 ENDDO
 ENDDO

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_tked(i,k,j) =0.25*(g_tke(i-1,k,j) +g_tke(i,k,j) +g_tke(i-1,k,j-1) &
 +g_tke(i,k,j-1))
 tked(i,k,j) =0.25*(tke(i-1,k,j) +tke(i,k,j) +tke(i-1,k,j-1) +tke(i,k,j-1))

 g_smnsmnd(i,k,j) =0.25*(g_smnsmn(i-1,k,j) +g_smnsmn(i,k,j) +g_smnsmn(i-1, &
 k,j-1) +g_smnsmn(i,k,j-1))
 smnsmnd(i,k,j) =0.25*(smnsmn(i-1,k,j) +smnsmn(i,k,j) +smnsmn(i-1,k,j-1) +smnsmn(i,k,j-1))

 g_ss11d(i,k,j) =0.25*(g_ss11(i-1,k,j) +g_ss11(i,k,j) +g_ss11(i-1,k,j-1) &
 +g_ss11(i,k,j-1))
 ss11d(i,k,j) =0.25*(ss11(i-1,k,j) +ss11(i,k,j) +ss11(i-1,k,j-1) +ss11(i,k,j-1))

 g_ss22d(i,k,j) =0.25*(g_ss22(i-1,k,j) +g_ss22(i,k,j) +g_ss22(i-1,k,j-1) &
 +g_ss22(i,k,j-1))
 ss22d(i,k,j) =0.25*(ss22(i-1,k,j) +ss22(i,k,j) +ss22(i-1,k,j-1) +ss22(i,k,j-1))

 g_ss13d(i,k,j) =0.25*(g_ss13(i,k+1,j) +g_ss13(i,k+1,j-1) +g_ss13(i,k,j) &
 +g_ss13(i,k,j-1))
 ss13d(i,k,j) =0.25*(ss13(i,k+1,j) +ss13(i,k+1,j-1) +ss13(i,k,j) +ss13(i,k,j-1))

 g_rr13d(i,k,j) =0.25*(g_rr13(i,k+1,j) +g_rr13(i,k+1,j-1) +g_rr13(i,k,j) &
 +g_rr13(i,k,j-1))
 rr13d(i,k,j) =0.25*(rr13(i,k+1,j) +rr13(i,k+1,j-1) +rr13(i,k,j) +rr13(i,k,j-1))

 g_ss23d(i,k,j) =0.25*(g_ss23(i,k+1,j) +g_ss23(i-1,k+1,j) +g_ss23(i,k,j) &
 +g_ss23(i-1,k,j))
 ss23d(i,k,j) =0.25*(ss23(i,k+1,j) +ss23(i-1,k+1,j) +ss23(i,k,j) +ss23(i-1,k,j))

 g_rr23d(i,k,j) =0.25*(g_rr23(i,k+1,j) +g_rr23(i-1,k+1,j) +g_rr23(i,k,j) &
 +g_rr23(i-1,k,j))
 rr23d(i,k,j) =0.25*(rr23(i,k+1,j) +rr23(i-1,k+1,j) +rr23(i,k,j) +rr23(i-1,k,j))

 ENDDO
 ENDDO
 ENDDO

 IF( config_flags%sfs_opt .EQ. 1 ) THEN

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0*2.0*(cs*g_delta)*(cs*delta)
 a =-1.0*(cs*delta)**2

 g_Tmpv1 =2.0*sqrt(2.0*smnsmnd(i,k,j))*g_ss12(i,k,j) +2.0*g_Sqrt(2.0* &
 g_smnsmnd(i,k,j), 2.0*smnsmnd(i,k,j))*ss12(i,k,j)
 Tmpv1 =2.0*sqrt(2.0*smnsmnd(i,k,j))*ss12(i,k,j)

 g_Tmpv2 =ss11d(i,k,j)*g_ss12(i,k,j) +g_ss11d(i,k,j)*ss12(i,k,j) 
 Tmpv2 =ss11d(i,k,j)*ss12(i,k,j)

 g_Tmpv3 =ss22d(i,k,j)*g_ss12(i,k,j) +g_ss22d(i,k,j)*ss12(i,k,j) 
 Tmpv3 =ss22d(i,k,j)*ss12(i,k,j)

 g_Tmpv4 =ss13d(i,k,j)*g_ss23d(i,k,j) +g_ss13d(i,k,j)*ss23d(i,k,j) 
 Tmpv4 =ss13d(i,k,j)*ss23d(i,k,j)

 g_Tmpv5 =ss11d(i,k,j)*g_rr12(i,k,j) +g_ss11d(i,k,j)*rr12(i,k,j) 
 Tmpv5 =ss11d(i,k,j)*rr12(i,k,j)

 g_Tmpv6 =ss13d(i,k,j)*g_rr23d(i,k,j) +g_ss13d(i,k,j)*rr23d(i,k,j) 
 Tmpv6 =ss13d(i,k,j)*rr23d(i,k,j)

 g_Tmpv7 =ss22d(i,k,j)*g_rr12(i,k,j) +g_ss22d(i,k,j)*rr12(i,k,j) 
 Tmpv7 =ss22d(i,k,j)*rr12(i,k,j)

 g_Tmpv8 =ss23d(i,k,j)*g_rr13d(i,k,j) +g_ss23d(i,k,j)*rr13d(i,k,j) 
 Tmpv8 =ss23d(i,k,j)*rr13d(i,k,j)

 g_Tmpv9 =a*(g_Tmpv1 +c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4) +c2*(g_Tmpv5 - &
 g_Tmpv6 -g_Tmpv7 -g_Tmpv8)) +g_a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4) &
 +c2*(Tmpv5 -Tmpv6 -Tmpv7 -Tmpv8)) 
 Tmpv9 =a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 -Tmpv6 -Tmpv7 -Tmpv8))

 g_m12(i,k,j) =g_Tmpv9
 m12(i,k,j) =Tmpv9

 ENDDO
 ENDDO
 ENDDO

 ELSE

 DO j =j_start,j_end
 DO k =kts,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0 *ce*g_delta
 a =-1.0 *ce*delta

 g_b =c3*g_delta
 b =c3*delta

 g_Tmpv1 =2.0*sqrt(tked(i,k,j))*g_s12(i,k,j) +2.0*g_Sqrt(g_tked(i,k,j) &
, tked(i,k,j))*s12(i,k,j) 
 Tmpv1 =2.0*sqrt(tked(i,k,j))*s12(i,k,j)

 g_Tmpv2 =ss11d(i,k,j)*g_ss12(i,k,j) +g_ss11d(i,k,j)*ss12(i,k,j) 
 Tmpv2 =ss11d(i,k,j)*ss12(i,k,j)

 g_Tmpv3 =ss22d(i,k,j)*g_ss12(i,k,j) +g_ss22d(i,k,j)*ss12(i,k,j) 
 Tmpv3 =ss22d(i,k,j)*ss12(i,k,j)

 g_Tmpv4 =ss13d(i,k,j)*g_ss23d(i,k,j) +g_ss13d(i,k,j)*ss23d(i,k,j) 
 Tmpv4 =ss13d(i,k,j)*ss23d(i,k,j)

 g_Tmpv5 =ss11d(i,k,j)*g_rr12(i,k,j) +g_ss11d(i,k,j)*rr12(i,k,j) 
 Tmpv5 =ss11d(i,k,j)*rr12(i,k,j)

 g_Tmpv6 =ss13d(i,k,j)*g_rr23d(i,k,j) +g_ss13d(i,k,j)*rr23d(i,k,j) 
 Tmpv6 =ss13d(i,k,j)*rr23d(i,k,j)

 g_Tmpv7 =ss22d(i,k,j)*g_rr12(i,k,j) +g_ss22d(i,k,j)*rr12(i,k,j) 
 Tmpv7 =ss22d(i,k,j)*rr12(i,k,j)

 g_Tmpv8 =ss23d(i,k,j)*g_rr13d(i,k,j) +g_ss23d(i,k,j)*rr13d(i,k,j) 
 Tmpv8 =ss23d(i,k,j)*rr13d(i,k,j)

 g_Tmpv9 =b*(c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4) +c2*(g_Tmpv5 -g_Tmpv6 - &
 g_Tmpv7 -g_Tmpv8)) +g_b*(c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 -Tmpv6 -Tmpv7 -Tmpv8)) 
 Tmpv9 =b*(c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 -Tmpv6 -Tmpv7 -Tmpv8))

 g_Tmpv10 =a*(g_Tmpv1 +g_Tmpv9) +g_a*(Tmpv1 +Tmpv9) 
 Tmpv10 =a*(Tmpv1 +Tmpv9)

 g_m12(i,k,j) =g_Tmpv10
 m12(i,k,j) =Tmpv10

 ENDDO
 ENDDO
 ENDDO
 ENDIF
 Return

 END SUBROUTINE g_calc_m12

 SUBROUTINE g_calc_m13(m13,g_m13,s11,g_s11,s33,g_s33,s12,g_s12,s13, &
 g_s13,s23,g_s23,r12,g_r12,r13,g_r13,r23,g_r23,smnsmn,g_smnsmn,tke, &
 g_tke,rdzw,g_rdzw,dx,dy,fnm,fnp,config_flags,ids,ide,jds,jde,kds,kde,ims,ime, &
 jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,its,ite,jts,jte,kts,kte)

 IMPLICIT NONE

 REAL :: Tmpv1,g_Tmpv1,Tmpv2,g_Tmpv2,Tmpv3,g_Tmpv3,Tmpv4,g_Tmpv4,Tmpv5, &
 g_Tmpv5,Tmpv6,g_Tmpv6,Tmpv7,g_Tmpv7,Tmpv8,g_Tmpv8,Tmpv9,g_Tmpv9, &
 Tmpv10,g_Tmpv10

   REAL :: g_Sqrt
 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: m13,g_m13

 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: s11,g_s11,s33,g_s33,s12,g_s12, &
 s13,g_s13,s23,g_s23,r12,g_r12,r13,g_r13,r23,g_r23,smnsmn,g_smnsmn, &
 tke,g_tke,rdzw,g_rdzw

 REAL :: dx,dy

 REAL,DIMENSION(kms:kme) :: fnm,fnp
 TYPE(grid_config_rec_type) :: config_flags
 INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe, &
 its,ite,jts,jte,kts,kte

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: ss11,g_ss11,ss33,g_ss33, &
 ss12,g_ss12,ss13,g_ss13,ss23,g_ss23,rr12,g_rr12,rr13,g_rr13,rr23,g_rr23

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: tkee,g_tkee,ss11e,g_ss11e, &
 ss33e,g_ss33e,ss12e,g_ss12e,ss23e,g_ss23e,rr12e,g_rr12e,rr23e,g_rr23e, &
 smnsmne,g_smnsmne
 REAL :: delta,g_delta,a,g_a,b,g_b
 INTEGER :: i,j,k,i_start,i_end,j_start,j_end,ktf,ie_ext

 ktf =min(kte,kde-1)

 i_start =its

 i_end =ite

 j_start =jts

 j_end =min(jte,jde-1)

 IF( config_flags%open_xs .OR. config_flags%specified .OR.   &
         config_flags%nested) i_start =max(ids+1,its)

 IF( config_flags%open_xe .OR. config_flags%specified .OR.   &
         config_flags%nested) i_end =min(ide-1,ite)

 IF( config_flags%open_ys .OR. config_flags%specified .OR.   &
         config_flags%nested) j_start =max(jds+1,jts)

 IF( config_flags%open_ye .OR. config_flags%specified .OR.   &
         config_flags%nested) j_end =min(jde-2,jte)

 IF( config_flags%periodic_x ) i_start =its

 IF( config_flags%periodic_x ) i_end =ite

 ie_ext =1

 i_end =i_end+ie_ext

 DO j =j_start,j_end+1
 DO k =kts,ktf
 DO i =i_start-1,i_end

 g_ss11(i,k,j) =g_s11(i,k,j)/2.0
 ss11(i,k,j) =s11(i,k,j)/2.0

 g_ss33(i,k,j) =g_s33(i,k,j)/2.0
 ss33(i,k,j) =s33(i,k,j)/2.0

 g_ss12(i,k,j) =g_s12(i,k,j)/2.0
 ss12(i,k,j) =s12(i,k,j)/2.0

 g_ss13(i,k,j) =g_s13(i,k,j)/2.0
 ss13(i,k,j) =s13(i,k,j)/2.0

 g_ss23(i,k,j) =g_s23(i,k,j)/2.0
 ss23(i,k,j) =s23(i,k,j)/2.0

 g_rr12(i,k,j) =g_r12(i,k,j)/2.0
 rr12(i,k,j) =r12(i,k,j)/2.0

 g_rr13(i,k,j) =g_r13(i,k,j)/2.0
 rr13(i,k,j) =r13(i,k,j)/2.0

 g_rr23(i,k,j) =g_r23(i,k,j)/2.0
 rr23(i,k,j) =r23(i,k,j)/2.0

 ENDDO
 ENDDO
 ENDDO

 DO j =j_start,j_end
 DO k =kts+1,ktf
 DO i =i_start,i_end

 g_tkee(i,k,j) =0.5*(fnm(k)*(g_tke(i,k,j) +g_tke(i-1,k,j)) +fnp(k) &
*(g_tke(i,k-1,j) +g_tke(i-1,k-1,j)))
 tkee(i,k,j) =0.5*(fnm(k)*(tke(i,k,j) +tke(i-1,k,j)) +fnp(k)*(tke(i,k-1,j) +tke(i-1,k-1,j)))

 g_smnsmne(i,k,j) =0.5*(fnm(k)*(g_smnsmn(i,k,j) +g_smnsmn(i-1,k,j)) +fnp(k) &
*(g_smnsmn(i,k-1,j) +g_smnsmn(i-1,k-1,j)))
 smnsmne(i,k,j) =0.5*(fnm(k)*(smnsmn(i,k,j) +smnsmn(i-1,k,j)) +fnp(k)*(smnsmn(i,k-1,j) &
 +smnsmn(i-1,k-1,j)))

 g_ss11e(i,k,j) =0.5*(fnm(k)*(g_ss11(i,k,j) +g_ss11(i-1,k,j)) +fnp(k) &
*(g_ss11(i,k-1,j) +g_ss11(i-1,k-1,j)))
 ss11e(i,k,j) =0.5*(fnm(k)*(ss11(i,k,j) +ss11(i-1,k,j)) +fnp(k)*(ss11(i,k-1,j) &
 +ss11(i-1,k-1,j)))

 g_ss33e(i,k,j) =0.5*(fnm(k)*(g_ss33(i,k,j) +g_ss33(i-1,k,j)) +fnp(k) &
*(g_ss33(i,k-1,j) +g_ss33(i-1,k-1,j)))
 ss33e(i,k,j) =0.5*(fnm(k)*(ss33(i,k,j) +ss33(i-1,k,j)) +fnp(k)*(ss33(i,k-1,j) &
 +ss33(i-1,k-1,j)))

 g_ss12e(i,k,j) =0.5*(fnm(k)*(g_ss12(i,k,j) +g_ss12(i,k,j+1)) +fnp(k) &
*(g_ss12(i,k-1,j) +g_ss12(i,k-1,j+1)))
 ss12e(i,k,j) =0.5*(fnm(k)*(ss12(i,k,j) +ss12(i,k,j+1)) +fnp(k)*(ss12(i,k-1,j) &
 +ss12(i,k-1,j+1)))

 g_rr12e(i,k,j) =0.5*(fnm(k)*(g_rr12(i,k,j) +g_rr12(i,k,j+1)) +fnp(k) &
*(g_rr12(i,k-1,j) +g_rr12(i,k-1,j+1)))
 rr12e(i,k,j) =0.5*(fnm(k)*(rr12(i,k,j) +rr12(i,k,j+1)) +fnp(k)*(rr12(i,k-1,j) &
 +rr12(i,k-1,j+1)))

 g_ss23e(i,k,j) =0.25*(g_ss23(i,k,j) +g_ss23(i,k,j+1) +g_ss23(i-1,k,j) &
 +g_ss23(i-1,k,j+1))
 ss23e(i,k,j) =0.25*(ss23(i,k,j) +ss23(i,k,j+1) +ss23(i-1,k,j) +ss23(i-1,k,j+1))

 g_rr23e(i,k,j) =0.25*(g_rr23(i,k,j) +g_rr23(i,k,j+1) +g_rr23(i-1,k,j) &
 +g_rr23(i-1,k,j+1))
 rr23e(i,k,j) =0.25*(rr23(i,k,j) +rr23(i,k,j+1) +rr23(i-1,k,j) +rr23(i-1,k,j+1))

 ENDDO
 ENDDO
 ENDDO

 IF( config_flags%sfs_opt .EQ. 1 ) THEN

 DO j =j_start,j_end
 DO k =kts+1,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0*2.0*(cs*g_delta)*(cs*delta)
 a =-1.0*(cs*delta)**2

 g_Tmpv1 =2.0*sqrt(2.0*smnsmne(i,k,j))*g_ss13(i,k,j) +2.0*g_Sqrt(2.0* &
 g_smnsmne(i,k,j), 2.0*smnsmne(i,k,j))*ss13(i,k,j) 
 Tmpv1 =2.0*sqrt(2.0*smnsmne(i,k,j))*ss13(i,k,j)

 g_Tmpv2 =ss11e(i,k,j)*g_ss13(i,k,j) +g_ss11e(i,k,j)*ss13(i,k,j) 
 Tmpv2 =ss11e(i,k,j)*ss13(i,k,j)

 g_Tmpv3 =ss12e(i,k,j)*g_ss23e(i,k,j) +g_ss12e(i,k,j)*ss23e(i,k,j) 
 Tmpv3 =ss12e(i,k,j)*ss23e(i,k,j)

 g_Tmpv4 =ss13(i,k,j)*g_ss33e(i,k,j) +g_ss13(i,k,j)*ss33e(i,k,j) 
 Tmpv4 =ss13(i,k,j)*ss33e(i,k,j)

 g_Tmpv5 =ss11e(i,k,j)*g_rr13(i,k,j) +g_ss11e(i,k,j)*rr13(i,k,j) 
 Tmpv5 =ss11e(i,k,j)*rr13(i,k,j)

 g_Tmpv6 =ss12e(i,k,j)*g_rr23e(i,k,j) +g_ss12e(i,k,j)*rr23e(i,k,j) 
 Tmpv6 =ss12e(i,k,j)*rr23e(i,k,j)

 g_Tmpv7 =ss23e(i,k,j)*g_rr12e(i,k,j) +g_ss23e(i,k,j)*rr12e(i,k,j) 
 Tmpv7 =ss23e(i,k,j)*rr12e(i,k,j)

 g_Tmpv8 =ss33e(i,k,j)*g_rr13(i,k,j) +g_ss33e(i,k,j)*rr13(i,k,j) 
 Tmpv8 =ss33e(i,k,j)*rr13(i,k,j)

 g_Tmpv9 =a*(g_Tmpv1 +c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4) +c2*(g_Tmpv5 + &
 g_Tmpv6 -g_Tmpv7 -g_Tmpv8)) +g_a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4) &
 +c2*(Tmpv5 +Tmpv6 -Tmpv7 -Tmpv8)) 
 Tmpv9 =a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 +Tmpv6 -Tmpv7 -Tmpv8))

 g_m13(i,k,j) =g_Tmpv9
 m13(i,k,j) =Tmpv9

 ENDDO
 ENDDO
 ENDDO

 ELSE

 DO j =j_start,j_end
 DO k =kts+1,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0 *ce*g_delta
 a =-1.0 *ce*delta

 g_b =c3*g_delta
 b =c3*delta

 g_Tmpv1 =2.0*sqrt(tkee(i,k,j))*g_ss13(i,k,j) +2.0*g_Sqrt(g_tkee(i,k,j) &
, tkee(i,k,j))*ss13(i,k,j) 
 Tmpv1 =2.0*sqrt(tkee(i,k,j))*ss13(i,k,j)

 g_Tmpv2 =ss11e(i,k,j)*g_ss13(i,k,j) +g_ss11e(i,k,j)*ss13(i,k,j) 
 Tmpv2 =ss11e(i,k,j)*ss13(i,k,j)

 g_Tmpv3 =ss12e(i,k,j)*g_ss23e(i,k,j) +g_ss12e(i,k,j)*ss23e(i,k,j) 
 Tmpv3 =ss12e(i,k,j)*ss23e(i,k,j)

 g_Tmpv4 =ss13(i,k,j)*g_ss33e(i,k,j) +g_ss13(i,k,j)*ss33e(i,k,j) 
 Tmpv4 =ss13(i,k,j)*ss33e(i,k,j)

 g_Tmpv5 =ss11e(i,k,j)*g_rr13(i,k,j) +g_ss11e(i,k,j)*rr13(i,k,j) 
 Tmpv5 =ss11e(i,k,j)*rr13(i,k,j)

 g_Tmpv6 =ss12e(i,k,j)*g_rr23e(i,k,j) +g_ss12e(i,k,j)*rr23e(i,k,j) 
 Tmpv6 =ss12e(i,k,j)*rr23e(i,k,j)

 g_Tmpv7 =ss23e(i,k,j)*g_rr12e(i,k,j) +g_ss23e(i,k,j)*rr12e(i,k,j) 
 Tmpv7 =ss23e(i,k,j)*rr12e(i,k,j)

 g_Tmpv8 =ss33e(i,k,j)*g_rr13(i,k,j) +g_ss33e(i,k,j)*rr13(i,k,j) 
 Tmpv8 =ss33e(i,k,j)*rr13(i,k,j)

 g_Tmpv9 =b*(c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4) +c2*(g_Tmpv5 +g_Tmpv6 - &
 g_Tmpv7 -g_Tmpv8)) +g_b*(c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 +Tmpv6 -Tmpv7 -Tmpv8)) 
 Tmpv9 =b*(c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 +Tmpv6 -Tmpv7 -Tmpv8))

 g_Tmpv10 =a*(g_Tmpv1 +g_Tmpv9) +g_a*(Tmpv1 +Tmpv9) 
 Tmpv10 =a*(Tmpv1 +Tmpv9)

 g_m13(i,k,j) =g_Tmpv10
 m13(i,k,j) =Tmpv10

 ENDDO
 ENDDO
 ENDDO
 ENDIF
 Return

 END SUBROUTINE g_calc_m13

 SUBROUTINE g_calc_m23(m23,g_m23,s22,g_s22,s33,g_s33,s12,g_s12,s13, &
 g_s13,s23,g_s23,r12,g_r12,r13,g_r13,r23,g_r23,smnsmn,g_smnsmn,tke, &
 g_tke,rdzw,g_rdzw,dx,dy,fnm,fnp,config_flags,ids,ide,jds,jde,kds,kde,ims,ime, &
 jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe,its,ite,jts,jte,kts,kte)

 IMPLICIT NONE

 REAL :: Tmpv1,g_Tmpv1,Tmpv2,g_Tmpv2,Tmpv3,g_Tmpv3,Tmpv4,g_Tmpv4,Tmpv5, &
 g_Tmpv5,Tmpv6,g_Tmpv6,Tmpv7,g_Tmpv7,Tmpv8,g_Tmpv8,Tmpv9,g_Tmpv9, &
 Tmpv10,g_Tmpv10

   REAL :: g_Sqrt
 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: m23,g_m23

 REAL,DIMENSION(ims:ime,kms:kme,jms:jme) :: s22,g_s22,s33,g_s33,s12,g_s12, &
 s13,g_s13,s23,g_s23,r12,g_r12,r13,g_r13,r23,g_r23,smnsmn,g_smnsmn, &
 tke,g_tke,rdzw,g_rdzw

 REAL :: dx,dy

 REAL,DIMENSION(kms:kme) :: fnm,fnp
 TYPE(grid_config_rec_type) :: config_flags
 INTEGER :: ids,ide,jds,jde,kds,kde,ims,ime,jms,jme,kms,kme,ips,ipe,jps,jpe,kps,kpe, &
 its,ite,jts,jte,kts,kte

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: ss22,g_ss22,ss33,g_ss33, &
 ss12,g_ss12,ss13,g_ss13,ss23,g_ss23,rr12,g_rr12,rr13,g_rr13,rr23,g_rr23

 REAL,DIMENSION(its-1:ite+1,kms:kme,jts-1:jte+1) :: tkef,g_tkef,ss22f,g_ss22f, &
 ss33f,g_ss33f,ss12f,g_ss12f,ss13f,g_ss13f,rr12f,g_rr12f,rr13f,g_rr13f, &
 smnsmnf,g_smnsmnf
 REAL :: delta,g_delta,a,g_a,b,g_b
 INTEGER :: i,j,k,i_start,i_end,j_start,j_end,ktf,je_ext

 ktf =min(kte,kde-1)

 i_start =its

 i_end =min(ite,ide-1)

 j_start =jts

 j_end =jte

 IF( config_flags%open_xs .OR. config_flags%specified .OR.   &
         config_flags%nested) i_start =max(ids+1,its)

 IF( config_flags%open_xe .OR. config_flags%specified .OR.   &
         config_flags%nested) i_end =min(ide-2,ite)

 IF( config_flags%open_ys .OR. config_flags%specified .OR.   &
         config_flags%nested) j_start =max(jds+1,jts)

 IF( config_flags%open_ye .OR. config_flags%specified .OR.   &
         config_flags%nested) j_end =min(jde-1,jte)

 IF( config_flags%periodic_x ) i_start =its

 IF( config_flags%periodic_x ) i_end =min(ite,ide-1)

 je_ext =1

 j_end =j_end+je_ext

 DO j =j_start-1,j_end
 DO k =kts,ktf
 DO i =i_start,i_end+1

 g_ss22(i,k,j) =g_s22(i,k,j)/2.0
 ss22(i,k,j) =s22(i,k,j)/2.0

 g_ss33(i,k,j) =g_s33(i,k,j)/2.0
 ss33(i,k,j) =s33(i,k,j)/2.0

 g_ss12(i,k,j) =g_s12(i,k,j)/2.0
 ss12(i,k,j) =s12(i,k,j)/2.0

 g_ss13(i,k,j) =g_s13(i,k,j)/2.0
 ss13(i,k,j) =s13(i,k,j)/2.0

 g_ss23(i,k,j) =g_s23(i,k,j)/2.0
 ss23(i,k,j) =s23(i,k,j)/2.0

 g_rr12(i,k,j) =g_r12(i,k,j)/2.0
 rr12(i,k,j) =r12(i,k,j)/2.0

 g_rr13(i,k,j) =g_r13(i,k,j)/2.0
 rr13(i,k,j) =r13(i,k,j)/2.0

 g_rr23(i,k,j) =g_r23(i,k,j)/2.0
 rr23(i,k,j) =r23(i,k,j)/2.0

 ENDDO
 ENDDO
 ENDDO

 DO j =j_start,j_end
 DO k =kts+1,ktf
 DO i =i_start,i_end

 g_tkef(i,k,j) =0.5*(fnm(k)*(g_tke(i,k,j) +g_tke(i,k,j-1)) +fnp(k) &
*(g_tke(i,k-1,j) +g_tke(i,k-1,j-1)))
 tkef(i,k,j) =0.5*(fnm(k)*(tke(i,k,j) +tke(i,k,j-1)) +fnp(k)*(tke(i,k-1,j) +tke(i,k-1,j-1)))

 g_smnsmnf(i,k,j) =0.5*(fnm(k)*(g_smnsmn(i,k,j) +g_smnsmn(i,k,j-1)) +fnp(k) &
*(g_smnsmn(i,k-1,j) +g_smnsmn(i,k-1,j-1)))
 smnsmnf(i,k,j) =0.5*(fnm(k)*(smnsmn(i,k,j) +smnsmn(i,k,j-1)) +fnp(k)*(smnsmn(i,k-1,j) &
 +smnsmn(i,k-1,j-1)))

 g_ss22f(i,k,j) =0.5*(fnm(k)*(g_ss22(i,k,j) +g_ss22(i,k,j-1)) +fnp(k) &
*(g_ss22(i,k-1,j) +g_ss22(i,k-1,j-1)))
 ss22f(i,k,j) =0.5*(fnm(k)*(ss22(i,k,j) +ss22(i,k,j-1)) +fnp(k)*(ss22(i,k-1,j) &
 +ss22(i,k-1,j-1)))

 g_ss33f(i,k,j) =0.5*(fnm(k)*(g_ss33(i,k,j) +g_ss33(i,k,j-1)) +fnp(k) &
*(g_ss33(i,k-1,j) +g_ss33(i,k-1,j-1)))
 ss33f(i,k,j) =0.5*(fnm(k)*(ss33(i,k,j) +ss33(i,k,j-1)) +fnp(k)*(ss33(i,k-1,j) &
 +ss33(i,k-1,j-1)))

 g_ss12f(i,k,j) =0.5*(fnm(k)*(g_ss12(i,k,j) +g_ss12(i+1,k,j)) +fnp(k) &
*(g_ss12(i,k-1,j) +g_ss12(i+1,k-1,j)))
 ss12f(i,k,j) =0.5*(fnm(k)*(ss12(i,k,j) +ss12(i+1,k,j)) +fnp(k)*(ss12(i,k-1,j) &
 +ss12(i+1,k-1,j)))

 g_rr12f(i,k,j) =0.5*(fnm(k)*(g_rr12(i,k,j) +g_rr12(i+1,k,j)) +fnp(k) &
*(g_rr12(i,k-1,j) +g_rr12(i+1,k-1,j)))
 rr12f(i,k,j) =0.5*(fnm(k)*(rr12(i,k,j) +rr12(i+1,k,j)) +fnp(k)*(rr12(i,k-1,j) &
 +rr12(i+1,k-1,j)))

 g_ss13f(i,k,j) =0.25*(g_ss13(i,k,j) +g_ss13(i,k,j-1) +g_ss13(i+1,k,j-1) &
 +g_ss13(i+1,k,j))
 ss13f(i,k,j) =0.25*(ss13(i,k,j) +ss13(i,k,j-1) +ss13(i+1,k,j-1) +ss13(i+1,k,j))

 g_rr13f(i,k,j) =0.25*(g_rr13(i,k,j) +g_rr13(i,k,j-1) +g_rr13(i+1,k,j-1) &
 +g_rr13(i+1,k,j))
 rr13f(i,k,j) =0.25*(rr13(i,k,j) +rr13(i,k,j-1) +rr13(i+1,k,j-1) +rr13(i+1,k,j))

 ENDDO
 ENDDO
 ENDDO

 IF( config_flags%sfs_opt .EQ. 1 ) THEN

 DO j =j_start,j_end
 DO k =kts+1,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0*2.0*(cs*g_delta)*(cs*delta)
 a =-1.0*(cs*delta)**2

 g_Tmpv1 =2.0*sqrt(2.0*smnsmnf(i,k,j))*g_ss23(i,k,j) +2.0*g_Sqrt(2.0* &
 g_smnsmnf(i,k,j), 2.0*smnsmnf(i,k,j))*ss23(i,k,j) 
 Tmpv1 =2.0*sqrt(2.0*smnsmnf(i,k,j))*ss23(i,k,j)

 g_Tmpv2 =ss12f(i,k,j)*g_ss13f(i,k,j) +g_ss12f(i,k,j)*ss13f(i,k,j) 
 Tmpv2 =ss12f(i,k,j)*ss13f(i,k,j)

 g_Tmpv3 =ss22f(i,k,j)*g_ss23(i,k,j) +g_ss22f(i,k,j)*ss23(i,k,j) 
 Tmpv3 =ss22f(i,k,j)*ss23(i,k,j)

 g_Tmpv4 =ss23(i,k,j)*g_ss33f(i,k,j) +g_ss23(i,k,j)*ss33f(i,k,j) 
 Tmpv4 =ss23(i,k,j)*ss33f(i,k,j)

 g_Tmpv5 =ss12f(i,k,j)*g_rr13f(i,k,j) +g_ss12f(i,k,j)*rr13f(i,k,j) 
 Tmpv5 =ss12f(i,k,j)*rr13f(i,k,j)

 g_Tmpv6 =ss22f(i,k,j)*g_rr23(i,k,j) +g_ss22f(i,k,j)*rr23(i,k,j) 
 Tmpv6 =ss22f(i,k,j)*rr23(i,k,j)

 g_Tmpv7 =ss13f(i,k,j)*g_rr12f(i,k,j) +g_ss13f(i,k,j)*rr12f(i,k,j) 
 Tmpv7 =ss13f(i,k,j)*rr12f(i,k,j)

 g_Tmpv8 =ss33f(i,k,j)*g_rr23(i,k,j) +g_ss33f(i,k,j)*rr23(i,k,j) 
 Tmpv8 =ss33f(i,k,j)*rr23(i,k,j)

 g_Tmpv9 =a*(g_Tmpv1 +c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4) +c2*(g_Tmpv5 + &
 g_Tmpv6 +g_Tmpv7 -g_Tmpv8)) +g_a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4) &
 +c2*(Tmpv5 +Tmpv6 +Tmpv7 -Tmpv8)) 
 Tmpv9 =a*(Tmpv1 +c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 +Tmpv6 +Tmpv7 -Tmpv8))

 g_m23(i,k,j) =g_Tmpv9
 m23(i,k,j) =Tmpv9

 ENDDO
 ENDDO
 ENDDO

 ELSE

 DO j =j_start,j_end
 DO k =kts+1,ktf
 DO i =i_start,i_end

 g_delta =0.33333333*(-dx *dy*g_rdzw(i,k,j)/(rdzw(i,k,j)*rdzw(i,k,j))) &
*(dx *dy/rdzw(i,k,j))**(0.33333333 -1.0)
 delta =(dx *dy/rdzw(i,k,j))**0.33333333

 g_a =-1.0 *ce*g_delta
 a =-1.0 *ce*delta

 g_b =c3*g_delta
 b =c3*delta

 g_Tmpv1 =2.0*sqrt(tkef(i,k,j))*g_ss23(i,k,j) +2.0*g_Sqrt(g_tkef(i,k,j) &
, tkef(i,k,j))*ss23(i,k,j) 
 Tmpv1 =2.0*sqrt(tkef(i,k,j))*ss23(i,k,j)

 g_Tmpv2 =ss12f(i,k,j)*g_ss13f(i,k,j) +g_ss12f(i,k,j)*ss13f(i,k,j) 
 Tmpv2 =ss12f(i,k,j)*ss13f(i,k,j)

 g_Tmpv3 =ss22f(i,k,j)*g_ss23(i,k,j) +g_ss22f(i,k,j)*ss23(i,k,j) 
 Tmpv3 =ss22f(i,k,j)*ss23(i,k,j)

 g_Tmpv4 =ss23(i,k,j)*g_ss33f(i,k,j) +g_ss23(i,k,j)*ss33f(i,k,j) 
 Tmpv4 =ss23(i,k,j)*ss33f(i,k,j)

 g_Tmpv5 =ss12f(i,k,j)*g_rr13f(i,k,j) +g_ss12f(i,k,j)*rr13f(i,k,j) 
 Tmpv5 =ss12f(i,k,j)*rr13f(i,k,j)

 g_Tmpv6 =ss22f(i,k,j)*g_rr23(i,k,j) +g_ss22f(i,k,j)*rr23(i,k,j) 
 Tmpv6 =ss22f(i,k,j)*rr23(i,k,j)

 g_Tmpv7 =ss13f(i,k,j)*g_rr12f(i,k,j) +g_ss13f(i,k,j)*rr12f(i,k,j) 
 Tmpv7 =ss13f(i,k,j)*rr12f(i,k,j)

 g_Tmpv8 =ss33f(i,k,j)*g_rr23(i,k,j) +g_ss33f(i,k,j)*rr23(i,k,j) 
 Tmpv8 =ss33f(i,k,j)*rr23(i,k,j)

 g_Tmpv9 =b*(c1*(g_Tmpv2 +g_Tmpv3 +g_Tmpv4) +c2*(g_Tmpv5 +g_Tmpv6 + &
 g_Tmpv7 -g_Tmpv8)) +g_b*(c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 +Tmpv6 +Tmpv7 -Tmpv8)) 
 Tmpv9 =b*(c1*(Tmpv2 +Tmpv3 +Tmpv4) +c2*(Tmpv5 +Tmpv6 +Tmpv7 -Tmpv8))

 g_Tmpv10 =a*(g_Tmpv1 +g_Tmpv9) +g_a*(Tmpv1 +Tmpv9) 
 Tmpv10 =a*(Tmpv1 +Tmpv9)

 g_m23(i,k,j) =g_Tmpv10
 m23(i,k,j) =Tmpv10

 ENDDO
 ENDDO
 ENDDO
 ENDIF
 Return

 END SUBROUTINE g_calc_m23

 END MODULE g_module_sfs_nba

! REAL Function g_Sqrt(g_x,x)
!
! REAL g_x,x
!
! IF(x.GT.0.0) THEN 
!   g_Sqrt =0.5*g_x/sqrt(x) 
! ELSE 
!   Print*,'' 
!   Print*,'g_Sqrt is incorrectly evaluated by 0!' 
!   Print*,'Aborted from calc_m23' 
!   g_Sqrt =0.0 
! END IF

! RETURN 
! END