Changeset 3764 for branches/2012/dev_MERGE_2012/NEMOGCM/NEMO/OPA_SRC/SBC/sbcdcy.F90

Timestamp:

2013-01-23T15:33:04+01:00 (11 years ago)

Author:

smasson

Message:

dev_MERGE_2012: report bugfixes done in the trunk from r3555 to r3763 into dev_MERGE_2012

File:

• branches/2012/dev_MERGE_2012/NEMOGCM/NEMO/OPA_SRC/SBC/sbcdcy.F90 (modified) (7 diffs)

branches/2012/dev_MERGE_2012/NEMOGCM/NEMO/OPA_SRC/SBC/sbcdcy.F90

@@ -105 +105 @@
          rcc(:,:) = zconvrad * glamt(:,:) - rpi
          ! time of midday
-         rtmd(:,:) = 0.5 - glamt(:,:) / 360.
-         rtmd(:,:) = MOD( (rtmd(:,:) + 1.), 1. )
+         rtmd(:,:) = 0.5_wp - glamt(:,:) / 360._wp
+         rtmd(:,:) = MOD( (rtmd(:,:) + 1._wp) , 1._wp)
       ENDIF
 
@@ -121 +121 @@
          zdsws = REAL(11 + nday_year, wp)
          ! declination of the earths orbit
-         zdecrad = (-23.5 * zconvrad) * COS( zdsws * ztwopi / REAL(nyear_len(1),wp) )
+         zdecrad = (-23.5_wp * zconvrad) * COS( zdsws * ztwopi / REAL(nyear_len(1),wp) )
          ! Compute A and B needed to compute the time integral of the diurnal cycle
 
@@ -138 +138 @@
          DO jj = 1, jpj
             DO ji = 1, jpi
-               IF ( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h
+               IF ( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h
          ! When is it night?
                   ztx = zinvtwopi * (ACOS(rab(ji,jj)) - rcc(ji,jj))
                   ztest = -rbb(ji,jj) * SIN( rcc(ji,jj) + ztwopi * ztx )
          ! is it dawn or dusk?
-                  IF ( ztest > 0 ) THEN
+                  IF ( ztest > 0._wp ) THEN
                      rdawn(ji,jj) = ztx
                      rdusk(ji,jj) = rtmd(ji,jj) + ( rtmd(ji,jj) - rdawn(ji,jj) )
@@ -151 +151 @@
                   ENDIF
                ELSE
-                  rdawn(ji,jj) = rtmd(ji,jj) + 0.5
+                  rdawn(ji,jj) = rtmd(ji,jj) + 0.5_wp
                   rdusk(ji,jj) = rdawn(ji,jj)
                ENDIF
@@ -158 +158 @@
          rdawn(:,:) = MOD( (rdawn(:,:) + 1._wp), 1._wp )
          rdusk(:,:) = MOD( (rdusk(:,:) + 1._wp), 1._wp )
-         !     2.2 Compute the scalling function:
-         !         S* = the inverse of the time integral of the diurnal cycle from dawm to dusk
+         !     2.2 Compute the scaling function:
+         !         S* = the inverse of the time integral of the diurnal cycle from dawn to dusk
+         !         Avoid possible infinite scaling factor, associated with very short daylight
+         !         periods, by ignoring periods less than 1/1000th of a day (ticket #1040)
          DO jj = 1, jpj
             DO ji = 1, jpi
-               IF ( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h
+               IF ( ABS(rab(ji,jj)) < 1._wp ) THEN         ! day duration is less than 24h
+                  rscal(ji,jj) = 0.0_wp
                   IF ( rdawn(ji,jj) < rdusk(ji,jj) ) THEN      ! day time in one part
-                     rscal(ji,jj) = fintegral(rdawn(ji,jj), rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
-                     rscal(ji,jj) = 1. / rscal(ji,jj)
+                     IF( (rdusk(ji,jj) - rdawn(ji,jj) ) .ge. 0.001_wp ) THEN
+                       rscal(ji,jj) = fintegral(rdawn(ji,jj), rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                       rscal(ji,jj) = 1._wp / rscal(ji,jj)
+                     ENDIF
                   ELSE                                         ! day time in two parts
-                     rscal(ji,jj) = fintegral(0., rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))   &
-                        &         + fintegral(rdawn(ji,jj), 1., raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
-                     rscal(ji,jj) = 1. / rscal(ji,jj)
+                     IF( (rdusk(ji,jj) + (1._wp - rdawn(ji,jj)) ) .ge. 0.001_wp ) THEN
+                       rscal(ji,jj) = fintegral(0._wp, rdusk(ji,jj), raa(ji,jj), rbb(ji,jj), rcc(ji,jj))   &
+                          &         + fintegral(rdawn(ji,jj), 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                       rscal(ji,jj) = 1. / rscal(ji,jj)
+                     ENDIF
                   ENDIF
                ELSE
                   IF ( raa(ji,jj) > rbb(ji,jj) ) THEN         ! 24h day
-                     rscal(ji,jj) = fintegral(0., 1., raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
-                     rscal(ji,jj) = 1. / rscal(ji,jj)
+                     rscal(ji,jj) = fintegral(0._wp, 1._wp, raa(ji,jj), rbb(ji,jj), rcc(ji,jj)) 
+                     rscal(ji,jj) = 1._wp / rscal(ji,jj)
                   ELSE                                          ! No day
-                     rscal(ji,jj) = 0.e0
+                     rscal(ji,jj) = 0.0_wp
                   ENDIF
                ENDIF
@@ -193 +200 @@
          DO ji = 1, jpi
             ztmpm = 0.0
-            IF( ABS(rab(ji,jj)) < 1 ) THEN         ! day duration is less than 24h
+            IF( ABS(rab(ji,jj)) < 1. ) THEN         ! day duration is less than 24h
                !
                IF( rdawn(ji,jj) < rdusk(ji,jj) ) THEN       ! day time in one part
@@ -227 +234 @@
                   !
                ELSE                                         ! No day
-                  zqsrout(ji,jj) = 0.e0
+                  zqsrout(ji,jj) = 0.0_wp
                   imask_night(ji,jj) = 1
                ENDIF