Changeset 5993 for branches/2015/dev_r5187_UKMO13_simplification/NEMOGCM/TOOLS/SCOORD_GEN/src/scoord_gen.F90
- Timestamp:
- 2015-12-03T14:55:36+01:00 (8 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5187_UKMO13_simplification/NEMOGCM/TOOLS/SCOORD_GEN/src/scoord_gen.F90
r5794 r5993 339 339 z_gsigw3 = 0. ; z_gsigt3 = 0. ; z_gsi3w3 = 0. 340 340 z_esigt3 = 0. ; z_esigw3 = 0. 341 z_esigt3p1 = 0. ; z_esigw3p1 = 0.342 341 z_esigtu3 = 0. ; z_esigtv3 = 0. ; z_esigtf3 = 0. 343 342 z_esigwu3 = 0. ; z_esigwv3 = 0. … … 664 663 665 664 END SUBROUTINE s_tanh 666 667 FUNCTION fssig( pk ) RESULT( pf )668 !!----------------------------------------------------------------------669 !! *** ROUTINE fssig ***670 !!671 !! ** Purpose : provide the analytical function in s-coordinate672 !!673 !! ** Method : the function provide the non-dimensional position of674 !! T and W (i.e. between 0 and 1)675 !! T-points at integer values (between 1 and jpk)676 !! W-points at integer values - 1/2 (between 0.5 and jpk-0.5)677 !!----------------------------------------------------------------------678 USE utils, ONLY : wp,rn_theta,rn_thetb,jpk679 IMPLICIT NONE680 REAL(wp), INTENT(in) :: pk ! continuous "k" coordinate681 REAL(wp) :: pf ! sigma value682 !!----------------------------------------------------------------------683 !684 pf = ( TANH( rn_theta * ( -(pk-0.5) / REAL(jpk-1) + rn_thetb ) ) &685 & - TANH( rn_thetb * rn_theta ) ) &686 & * ( COSH( rn_theta ) &687 & + COSH( rn_theta * ( 2. * rn_thetb - 1. ) ) ) &688 & / ( 2. * SINH( rn_theta ) )689 !690 END FUNCTION fssig691 692 693 FUNCTION fssig1( pk1, pbb ) RESULT( pf1 )694 !!----------------------------------------------------------------------695 !! *** ROUTINE fssig1 ***696 !!697 !! ** Purpose : provide the Song and Haidvogel version of the analytical function in s-coordinate698 !!699 !! ** Method : the function provides the non-dimensional position of700 !! T and W (i.e. between 0 and 1)701 !! T-points at integer values (between 1 and jpk)702 !! W-points at integer values - 1/2 (between 0.5 and jpk-0.5)703 !!----------------------------------------------------------------------704 USE utils, ONLY : wp, jpk, rn_theta705 IMPLICIT NONE706 REAL(wp), INTENT(in) :: pk1 ! continuous "k" coordinate707 REAL(wp), INTENT(in) :: pbb ! Stretching coefficient708 REAL(wp) :: pf1 ! sigma value709 !!----------------------------------------------------------------------710 !711 IF ( rn_theta == 0 ) then ! uniform sigma712 pf1 = - ( pk1 - 0.5 ) / REAL( jpk-1 )713 ELSE ! stretched sigma714 pf1 = ( 1. - pbb ) * ( SINH( rn_theta*(-(pk1-0.5)/REAL(jpk-1)) ) ) / SINH( rn_theta ) &715 & + pbb * ( (TANH( rn_theta*( (-(pk1-0.5)/REAL(jpk-1)) + 0.5) ) - TANH( 0.5 * rn_theta ) ) &716 & / ( 2. * TANH( 0.5 * rn_theta ) ) )717 ENDIF718 !719 END FUNCTION fssig1720 721 722 FUNCTION fgamma( pk1, pzb, pzs, psmth) RESULT( p_gamma )723 !!----------------------------------------------------------------------724 !! *** ROUTINE fgamma ***725 !!726 !! ** Purpose : provide analytical function for the s-coordinate727 !!728 !! ** Method : the function provides the non-dimensional position of729 !! T and W (i.e. between 0 and 1)730 !! T-points at integer values (between 1 and jpk)731 !! W-points at integer values - 1/2 (between 0.5 and jpk-0.5)732 !!733 !! This method allows the maintenance of fixed surface and or734 !! bottom cell resolutions (cf. geopotential coordinates)735 !! within an analytically derived stretched S-coordinate framework.736 !!737 !! Reference : Siddorn and Furner, in prep738 !!----------------------------------------------------------------------739 USE utils, ONLY : jpk,wp,rn_alpha740 IMPLICIT NONE741 REAL(wp), INTENT(in ) :: pk1 ! continuous "k" coordinate742 REAL(wp) :: p_gamma ! stretched coordinate743 REAL(wp), INTENT(in ) :: pzb ! Bottom box depth744 REAL(wp), INTENT(in ) :: pzs ! surface box depth745 REAL(wp), INTENT(in ) :: psmth ! Smoothing parameter746 REAL(wp) :: za1,za2,za3 ! local variables747 REAL(wp) :: zn1,zn2 ! local variables748 REAL(wp) :: za,zb,zx ! local variables749 !!----------------------------------------------------------------------750 !751 752 zn1 = 1./(jpk-1.)753 zn2 = 1. - zn1754 755 za1 = (rn_alpha+2.0)*zn1**(rn_alpha+1.0)-(rn_alpha+1.0)*zn1**(rn_alpha+2.0)756 za2 = (rn_alpha+2.0)*zn2**(rn_alpha+1.0)-(rn_alpha+1.0)*zn2**(rn_alpha+2.0)757 za3 = (zn2**3.0 - za2)/( zn1**3.0 - za1)758 759 za = pzb - za3*(pzs-za1)-za2760 za = za/( zn2-0.5*(za2+zn2**2.0) - za3*(zn1-0.5*(za1+zn1**2.0) ) )761 zb = (pzs - za1 - za*( zn1-0.5*(za1+zn1**2.0 ) ) ) / (zn1**3.0 - za1)762 zx = 1.0-za/2.0-zb763 764 p_gamma = za*(pk1*(1.0-pk1/2.0))+zb*pk1**3.0 + &765 & zx*( (rn_alpha+2.0)*pk1**(rn_alpha+1.0)- &766 & (rn_alpha+1.0)*pk1**(rn_alpha+2.0) )767 p_gamma = p_gamma*psmth+pk1*(1.0-psmth)768 769 !770 END FUNCTION fgamma771
Note: See TracChangeset
for help on using the changeset viewer.