Changeset 12045


Ignore:
Timestamp:
2019-12-03T18:24:37+01:00 (8 weeks ago)
Author:
smueller
Message:

Modifications related to coding style and conventions, and addition of a reference in module tide_mod (ticket #2194)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/TDE/tide_mod.F90

    r12042 r12045  
    478478      sh_p = MOD( sh_p,  2*rpi ) 
    479479 
    480       cosI = 0.913694997 -0.035692561 *cos(sh_N) 
     480      cosI = 0.913694997_wp - 0.035692561_wp * COS( sh_N ) 
    481481 
    482482      sh_I = ACOS( cosI ) 
     
    485485      sh_tgn2 = tan(sh_N/2.0) 
    486486 
    487       at1=atan(1.01883*sh_tgn2) 
    488       at2=atan(0.64412*sh_tgn2) 
    489  
    490       sh_xi=-at1-at2+sh_N 
    491  
    492       IF( sh_N > rpi )   sh_xi=sh_xi-2.0*rpi 
     487      at1 = ATAN( 1.01883_wp * sh_tgn2 ) 
     488      at2 = ATAN( 0.64412_wp * sh_tgn2 ) 
     489 
     490      sh_xi = sh_N - at1 - at2 
     491 
     492      IF( sh_N > rpi ) sh_xi = sh_xi - 2.0_wp * rpi 
    493493 
    494494      sh_nu = at1 - at2 
    495495 
    496       !---------------------------------------------------------------------- 
    497       ! For constituents l2 k1 k2 
    498       !---------------------------------------------------------------------- 
    499  
    500       tgI2 = tan(sh_I/2.0) 
    501       P1   = sh_p-sh_xi 
    502  
    503       t2 = tgI2*tgI2 
    504       t4 = t2*t2 
    505       sh_x1ra = sqrt( 1.0-12.0*t2*cos(2.0*P1)+36.0*t4 ) 
    506  
    507       p = sin(2.0*P1) 
    508       q = 1.0/(6.0*t2)-cos(2.0*P1) 
    509       sh_R = atan(p/q) 
    510  
    511       p = sin(2.0*sh_I)*sin(sh_nu) 
    512       q = sin(2.0*sh_I)*cos(sh_nu)+0.3347 
    513       sh_nuprim = atan(p/q) 
    514  
    515       s2 = sin(sh_I)*sin(sh_I) 
    516       p  = s2*sin(2.0*sh_nu) 
    517       q  = s2*cos(2.0*sh_nu)+0.0727 
    518       sh_nusec = 0.5*atan(p/q) 
     496      ! For computation of tidal constituents L2 K1 K2 
     497      tgI2 = tan( sh_I / 2.0_wp ) 
     498      P1   = sh_p - sh_xi 
     499      ! 
     500      t2   = tgI2 * tgI2 
     501      t4   = t2 * t2 
     502      sh_x1ra = SQRT( 1.0 - 12.0 * t2 * COS( 2.0 * P1 ) + 36.0_wp * t4 ) 
     503      ! 
     504      p = SIN( 2.0_wp * P1 ) 
     505      q = 1.0_wp / ( 6.0_wp * t2 ) - COS( 2.0_wp * P1 ) 
     506      sh_R = ATAN( p / q ) 
     507      ! 
     508      p = SIN( 2.0_wp * sh_I ) * SIN( sh_nu ) 
     509      q = SIN( 2.0_wp * sh_I ) * COS( sh_nu ) + 0.3347_wp 
     510      sh_nuprim = ATAN( p / q ) 
     511      ! 
     512      s2 = SIN( sh_I ) * SIN( sh_I ) 
     513      p  = s2 * SIN( 2.0_wp * sh_nu ) 
     514      q  = s2 * COS( 2.0_wp * sh_nu ) + 0.0727_wp 
     515      sh_nusec = 0.5_wp * ATAN( p / q ) 
    519516      ! 
    520517   END SUBROUTINE astronomic_angle 
     
    595592   RECURSIVE FUNCTION nodal_factort( kformula ) RESULT( zf ) 
    596593      !!---------------------------------------------------------------------- 
    597       !!---------------------------------------------------------------------- 
    598       INTEGER, INTENT(in) :: kformula 
    599       ! 
    600       REAL(wp)            :: zf 
    601       REAL(wp)            :: zs, zf1, zf2 
    602       CHARACTER(LEN=3)    :: clformula 
     594      !!                  ***  FUNCTION nodal_factort  *** 
     595      !! 
     596      !! ** Purpose : Compute amplitude correction factors due to nodal motion 
     597      !! 
     598      !! ** Reference : 
     599      !!       S58) Schureman, P. (1958): Manual of Harmonic Analysis and 
     600      !!            Prediction of Tides (Revised (1940) Edition (Reprinted 1958 
     601      !!            with corrections). Reprinted June 2001). U.S. Department of 
     602      !!            Commerce, Coast and Geodetic Survey Special Publication 
     603      !!            No. 98. Washington DC, United States Government Printing 
     604      !!            Office. 317 pp. DOI: 10.25607/OBP-155. 
     605      !!---------------------------------------------------------------------- 
     606      INTEGER, INTENT(in) ::   kformula 
     607      ! 
     608      REAL(wp)            ::   zf 
     609      REAL(wp)            ::   zs, zf1, zf2 
     610      CHARACTER(LEN=3)    ::   clformula 
    603611      !!---------------------------------------------------------------------- 
    604612      ! 
    605613      SELECT CASE( kformula ) 
    606614      ! 
    607       CASE( 0 )                  !==  formule 0, solar waves 
     615      CASE( 0 )                  ! Formula 0, solar waves 
    608616         zf = 1.0 
    609617         ! 
    610       CASE( 1 )                  !==  formule 1, compound waves (78 x 78) 
    611          zf=nodal_factort(78) 
     618      CASE( 1 )                  ! Formula 1, compound waves (78 x 78) 
     619         zf=nodal_factort( 78 ) 
    612620         zf = zf * zf 
    613621         ! 
    614       CASE ( 4 )                 !==  formule 4,  compound waves (78 x 235)  
    615          zf1 = nodal_factort( 78) 
     622      CASE ( 4 )                 ! Formula 4,  compound waves (78 x 235)  
     623         zf1 = nodal_factort( 78 ) 
    616624         zf  = nodal_factort(235) 
    617625         zf  = zf1 * zf 
    618626         ! 
    619       CASE( 18 )                 !==  formule 18,  compound waves (78 x 78 x 78 ) 
    620          zf1 = nodal_factort(78) 
     627      CASE( 18 )                 ! Formula 18,  compound waves (78 x 78 x 78 ) 
     628         zf1 = nodal_factort( 78 ) 
    621629         zf  = zf1 * zf1 * zf1 
    622630         ! 
    623       CASE( 20 )                 !==  formula 20, compound waves ( 78 x 78 x 78 x 78 ) 
    624          zf1 = nodal_factort(78) 
     631      CASE( 20 )                 ! Formula 20, compound waves ( 78 x 78 x 78 x 78 ) 
     632         zf1 = nodal_factort( 78 ) 
    625633         zf  = zf1 * zf1 * zf1 * zf1 
    626634         ! 
    627       CASE( 73 )                 !==  formule 73 
    628          zs = sin(sh_I) 
    629          zf = (2./3.-zs*zs)/0.5021 
    630          ! 
    631       CASE( 74 )                 !==  formule 74 
    632          zs = sin(sh_I) 
    633          zf = zs * zs / 0.1578 
    634          ! 
    635       CASE( 75 )                 !==  formule 75 
    636          zs = cos(sh_I/2) 
    637          zf = sin(sh_I) * zs * zs / 0.3800 
    638          ! 
    639       CASE( 76 )                 !==  formule 76 
    640          zf = sin(2*sh_I) / 0.7214 
    641          ! 
    642       CASE( 78 )                 !==  formule 78 
    643          zs = cos(sh_I/2) 
    644          zf = zs * zs * zs * zs / 0.9154 
    645          ! 
    646       CASE( 149 )                !==  formule 149 
    647          zs = cos(sh_I/2) 
    648          zf = zs*zs*zs*zs*zs*zs / 0.8758 
    649          ! 
    650       CASE( 215 )                !==  formule 215 
    651          zs = cos(sh_I/2) 
    652          zf = zs*zs*zs*zs / 0.9154 * sh_x1ra 
    653          ! 
    654       CASE( 227 )                !==  formule 227  
    655          zs = sin(2*sh_I) 
    656          zf = sqrt( 0.8965*zs*zs+0.6001*zs*cos (sh_nu)+0.1006 ) 
    657          ! 
    658       CASE ( 235 )               !==  formule 235  
    659          zs = sin(sh_I) 
    660          zf = sqrt( 19.0444*zs*zs*zs*zs + 2.7702*zs*zs*cos(2*sh_nu) + .0981 ) 
     635      CASE( 73 )                 ! Formula 73 of S58 
     636         zs = SIN( sh_I ) 
     637         zf = ( 2.0_wp / 3.0_wp - zs * zs ) / 0.5021_wp 
     638         ! 
     639      CASE( 74 )                 ! Formula 74 of S58 
     640         zs = SIN(sh_I) 
     641         zf = zs * zs / 0.1578_wp 
     642         ! 
     643      CASE( 75 )                 ! Formula 75 of S58 
     644         zs = COS( sh_I / 2.0_wp ) 
     645         zf = SIN( sh_I ) * zs * zs / 0.3800_wp 
     646         ! 
     647      CASE( 76 )                 ! Formula 76 of S58 
     648         zf = SIN( 2.0_wp * sh_I ) / 0.7214_wp 
     649         ! 
     650      CASE( 78 )                 ! Formula 78 of S58 
     651         zs = COS( sh_I/2 ) 
     652         zf = zs * zs * zs * zs / 0.9154_wp 
     653         ! 
     654      CASE( 149 )                ! Formula 149 of S58 
     655         zs = COS( sh_I/2 ) 
     656         zf = zs * zs * zs * zs * zs * zs / 0.8758_wp 
     657         ! 
     658      CASE( 215 )                ! Formula 215 of S58 with typo correction (0.9154 instead of 0.9145) 
     659         zs = COS( sh_I/2 ) 
     660         zf = zs * zs * zs * zs / 0.9154_wp * sh_x1ra 
     661         ! 
     662      CASE( 227 )                ! Formula 227 of S58 
     663         zs = SIN( 2.0_wp * sh_I ) 
     664         zf = SQRT( 0.8965_wp * zs * zs + 0.6001_wp * zs * COS( sh_nu ) + 0.1006_wp ) 
     665         ! 
     666      CASE ( 235 )               ! Formula 235 of S58 
     667         zs = SIN( sh_I ) 
     668         zf = SQRT( 19.0444_wp * zs * zs * zs * zs + 2.7702_wp * zs * zs * cos( 2.0_wp * sh_nu ) + 0.0981_wp ) 
    661669         ! 
    662670      CASE DEFAULT 
     
    670678   FUNCTION dayjul( kyr, kmonth, kday ) 
    671679      !!---------------------------------------------------------------------- 
    672       !!  *** THIS ROUTINE COMPUTES THE JULIAN DAY (AS A REAL VARIABLE) 
    673       !!---------------------------------------------------------------------- 
    674       INTEGER,INTENT(in) ::   kyr, kmonth, kday 
    675       ! 
    676       INTEGER,DIMENSION(12) ::  idayt, idays 
    677       INTEGER  ::   inc, ji 
    678       REAL(wp) ::   dayjul, zyq 
    679       ! 
    680       DATA idayt/0.,31.,59.,90.,120.,151.,181.,212.,243.,273.,304.,334./ 
    681       !!---------------------------------------------------------------------- 
    682       ! 
    683       idays(1) = 0. 
    684       idays(2) = 31. 
    685       inc = 0. 
    686       zyq = MOD( kyr-1900. , 4. ) 
    687       IF( zyq == 0.)   inc = 1. 
     680      !!                   ***  FUNCTION dayjul  *** 
     681      !! 
     682      !! Purpose :  compute the Julian day 
     683      !!---------------------------------------------------------------------- 
     684      INTEGER,INTENT(in)    ::   kyr, kmonth, kday 
     685      ! 
     686      INTEGER,DIMENSION(12) ::   idayt = (/ 0, 31, 59, 90, 120, 151, 181, 212, 243, 273, 304, 334 /) 
     687      INTEGER,DIMENSION(12) ::   idays 
     688      INTEGER               ::   inc, ji, zyq 
     689      REAL(wp)              ::   dayjul 
     690      !!---------------------------------------------------------------------- 
     691      ! 
     692      idays(1) = 0 
     693      idays(2) = 31 
     694      inc = 0.0_wp 
     695      zyq = MOD( kyr - 1900 , 4 ) 
     696      IF( zyq == 0 ) inc = 1 
    688697      DO ji = 3, 12 
    689          idays(ji)=idayt(ji)+inc 
     698         idays(ji) = idayt(ji) + inc 
    690699      END DO 
    691       dayjul = idays(kmonth) + kday 
     700      dayjul = REAL( idays(kmonth) + kday, KIND=wp ) 
    692701      ! 
    693702   END FUNCTION dayjul 
Note: See TracChangeset for help on using the changeset viewer.