Changeset 12042


Ignore:
Timestamp:
2019-12-03T16:41:39+01:00 (8 weeks ago)
Author:
smueller
Message:

Reduction of the number of real-valued literals used in subroutine astronomic_angle of 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

    r11864 r12042  
    6767   REAL(wp) ::   sh_xi, sh_nu, sh_nuprim, sh_nusec, sh_R   ! 
    6868   REAL(wp) ::   sh_I, sh_x1ra, sh_N                       ! 
     69 
     70   ! Longitudes on 1 Jan 1900, 00h and angular velocities (units of deg and 
     71   ! deg/h, respectively. The values of these module variables have been copied 
     72   ! from subroutine astronomic_angle of the version of this module used in 
     73   ! version 4.0 of NEMO 
     74   REAL(wp) ::   rlon00_N  =  259.1560564_wp               ! Longitude of ascending lunar node 
     75   REAL(wp) ::   romega_N  = -.0022064139_wp 
     76   REAL(wp) ::   rlon00_T  =  180.0_wp                     ! Mean solar angle (GMT) 
     77   REAL(wp) ::   romega_T  =  15.0_wp 
     78   REAL(wp) ::   rlon00_h  =  280.1895014_wp               ! Mean solar Longitude 
     79   REAL(wp) ::   romega_h  =  .0410686387_wp 
     80   REAL(wp) ::   rlon00_s  =  277.0256206_wp               ! Mean lunar Longitude 
     81   REAL(wp) ::   romega_s  =  .549016532_wp 
     82   REAL(wp) ::   rlon00_p1 =  281.2208569_wp               ! Longitude of solar perigee 
     83   REAL(wp) ::   romega_p1 =  .000001961_wp 
     84   REAL(wp) ::   rlon00_p  =  334.3837214_wp               ! Longitude of lunar perigee 
     85   REAL(wp) ::   romega_p  =  .004641834_wp 
    6986 
    7087   !!---------------------------------------------------------------------- 
     
    423440   SUBROUTINE astronomic_angle(psec_day) 
    424441      !!---------------------------------------------------------------------- 
    425       !!  tj is time elapsed since 1st January 1900, 0 hour, counted in julian 
    426       !!  century (e.g. time in days divide by 36525) 
     442      !!                 ***  ROUTINE astronomic_angle  *** 
     443      !!                       
     444      !! ** Purpose : Compute astronomic angles 
    427445      !!---------------------------------------------------------------------- 
    428446      INTEGER  ::   psec_day !   Number of seconds from midnight 
    429447      REAL(wp) ::   cosI, p, q, t2, t4, sin2I, s2, tgI2, P1, sh_tgn2, at1, at2 
    430       REAL(wp) ::   zqy , zsy, zday, zdj, zhfrac 
    431       !!---------------------------------------------------------------------- 
    432       ! 
    433       zqy = AINT( (nyear-1901.)/4. ) 
    434       zsy = nyear - 1900. 
     448      REAL(wp) ::   zqy , zsy, zday, zdj, zhfrac, zt 
     449      !!---------------------------------------------------------------------- 
     450      ! 
     451      ! Computation of the time from 1 Jan 1900, 00h in years 
     452      zqy = AINT( (nyear - 1901.0_wp) / 4.0_wp ) 
     453      zsy = nyear - 1900.0_wp 
    435454      ! 
    436455      zdj  = dayjul( nyear, nmonth, nday ) 
    437       zday = zdj + zqy - 1. 
    438       ! 
    439       zhfrac = psec_day / 3600. 
    440       ! 
    441       !---------------------------------------------------------------------- 
    442       !  Sh_n Longitude of ascending lunar node 
    443       !---------------------------------------------------------------------- 
    444       sh_N=(259.1560564-19.328185764*zsy-.0529539336*zday-.0022064139*zhfrac)*rad 
    445       !---------------------------------------------------------------------- 
    446       ! T mean solar angle (Greenwhich time) 
    447       !---------------------------------------------------------------------- 
    448       sh_T=(180.+zhfrac*(360./24.))*rad 
    449       !---------------------------------------------------------------------- 
    450       ! h mean solar Longitude 
    451       !---------------------------------------------------------------------- 
    452       sh_h=(280.1895014-.238724988*zsy+.9856473288*zday+.0410686387*zhfrac)*rad 
    453       !---------------------------------------------------------------------- 
    454       ! s mean lunar Longitude 
    455       !---------------------------------------------------------------------- 
    456       sh_s=(277.0256206+129.38482032*zsy+13.176396768*zday+.549016532*zhfrac)*rad 
    457       !---------------------------------------------------------------------- 
    458       ! p1 Longitude of solar perigee 
    459       !---------------------------------------------------------------------- 
    460       sh_p1=(281.2208569+.01717836*zsy+.000047064*zday+.000001961*zhfrac)*rad 
    461       !---------------------------------------------------------------------- 
    462       ! p Longitude of lunar perigee 
    463       !---------------------------------------------------------------------- 
    464       sh_p=(334.3837214+40.66246584*zsy+.111404016*zday+.004641834*zhfrac)*rad 
    465  
    466       sh_N = MOD( sh_N ,2*rpi ) 
    467       sh_s = MOD( sh_s ,2*rpi ) 
    468       sh_h = MOD( sh_h, 2*rpi ) 
    469       sh_p = MOD( sh_p, 2*rpi ) 
    470       sh_p1= MOD( sh_p1,2*rpi ) 
     456      zday = zdj + zqy - 1.0_wp 
     457      ! 
     458      zhfrac = psec_day / 3600.0_wp 
     459      ! 
     460      zt = zsy * 365.0_wp * 24.0_wp + zday * 24.0_wp + zhfrac 
     461      ! 
     462      ! Longitude of ascending lunar node 
     463      sh_N  = ( rlon00_N  + romega_N * zt     ) * rad 
     464      sh_N = MOD( sh_N,  2*rpi ) 
     465      ! Mean solar angle (Greenwhich time) 
     466      sh_T  = ( rlon00_T  + romega_T * zhfrac ) * rad 
     467      ! Mean solar Longitude 
     468      sh_h  = ( rlon00_h  + romega_h * zt     ) * rad 
     469      sh_h = MOD( sh_h,  2*rpi ) 
     470      ! Mean lunar Longitude 
     471      sh_s  = ( rlon00_s  + romega_s * zt     ) * rad 
     472      sh_s = MOD( sh_s,  2*rpi ) 
     473      ! Longitude of solar perigee 
     474      sh_p1 = ( rlon00_p1 + romega_p1 * zt    ) * rad 
     475      sh_p1= MOD( sh_p1, 2*rpi ) 
     476      ! Longitude of lunar perigee 
     477      sh_p  = ( rlon00_p  + romega_p  * zt    ) * rad 
     478      sh_p = MOD( sh_p,  2*rpi ) 
    471479 
    472480      cosI = 0.913694997 -0.035692561 *cos(sh_N) 
Note: See TracChangeset for help on using the changeset viewer.