Changeset 2661 for LMDZ5/trunk/libf/phylmd/physiq_mod.F90

Timestamp:

Oct 11, 2016, 5:28:09 PM (8 years ago)

Author:

Laurent Fairhead

Message:

Bug correction on the ozone reading:
When read_climoz = 1 or 2, values for ozone are read from a file with 360 daily values
so that a special calendar is used to read in the proper values when using a 365/366
day calendar in the model. Problem was the same special calendar was used to define
when the ozone was supposed to be updated resulting in an update that became
desynchronised with the actual change of days in the simulation. (is this making sense?)
Took advantage of this correction to rewrite the ozone update mechanism
LF

File:

• LMDZ5/trunk/libf/phylmd/physiq_mod.F90 (modified) (2 diffs)

LMDZ5/trunk/libf/phylmd/physiq_mod.F90

@@ -1016 +1016 @@
     ! ascending order
 
-    integer, save:: co3i = 0
-    !     time index in NetCDF file of current ozone fields
-    !$OMP THREADPRIVATE(co3i) 
-
     integer ro3i
     !     required time index in NetCDF file for the ozone fields, between 1
@@ -1913 +1909 @@
     ! Prescrire l'ozone et calculer l'albedo sur l'ocean.
     !
-    if (read_climoz >= 1) then
-       ! Ozone from a file
-       ! Update required ozone index:
-       ro3i = int((days_elapsed + jh_cur - jh_1jan) / year_len * 360.) + 1
-       if (ro3i == 361) ro3i = 360
-       ! (This should never occur, except perhaps because of roundup
-       ! error. See documentation.)
-       if (ro3i /= co3i) then
-          ! Update ozone field:
-          if (read_climoz == 1) then
-             call regr_pr_av(ncid_climoz, (/"tro3"/), julien=ro3i, &
-                  press_in_edg=press_climoz, paprs=paprs, v3=wo)
-          else
-             ! read_climoz == 2
-             call regr_pr_av(ncid_climoz, (/"tro3         ", &
-                  "tro3_daylight"/), julien=ro3i, press_in_edg=press_climoz, &
-                  paprs=paprs, v3=wo)
-          end if
-          ! Convert from mole fraction of ozone to column density of ozone in a
-          ! cell, in kDU:
-          forall (l = 1: read_climoz) wo(:, :, l) = wo(:, :, l) * rmo3 / rmd &
-               * zmasse / dobson_u / 1e3
-          ! (By regridding ozone values for LMDZ only once every 360th of
-          ! year, we have already neglected the variation of pressure in one
-          ! 360th of year. So do not recompute "wo" at each time step even if
-          ! "zmasse" changes a little.)
-          co3i = ro3i
-       end if
-    ELSEIF (MOD(itap-1,lmt_pas) == 0) THEN
+    ! Update ozone if day change
+    IF (MOD(itap-1,lmt_pas) == 0) THEN
+      IF (read_climoz == 0) THEN
        ! Once per day, update ozone from Royer:
-
-       IF (solarlong0<-999.) then
-          ! Generic case with evolvoing season
-          zzz=real(days_elapsed+1)
-       ELSE IF (abs(solarlong0-1000.)<1.e-4) then
-          ! Particular case with annual mean insolation
-          zzz=real(90) ! could be revisited
-!          IF (read_climoz/=-1) THEN
-!             abort_message ='read_climoz=-1 is recommended when ' &
-!                  // 'solarlong0=1000.'
-!             CALL abort_physic (modname,abort_message,1)
-!          ENDIF
-       ELSE
+        IF (solarlong0<-999.) then
+           ! Generic case with evolvoing season
+           zzz=real(days_elapsed+1)
+        ELSE IF (abs(solarlong0-1000.)<1.e-4) then
+           ! Particular case with annual mean insolation
+           zzz=real(90) ! could be revisited
+           IF (read_climoz/=-1) THEN
+              abort_message ='read_climoz=-1 is recommended when ' &
+                   // 'solarlong0=1000.'
+              CALL abort_physic (modname,abort_message,1)
+           ENDIF
+        ELSE
           ! Case where the season is imposed with solarlong0
           zzz=real(90) ! could be revisited
-       ENDIF
-       wo(:,:,1)=ozonecm(latitude_deg, paprs,read_climoz,rjour=zzz)
+        ENDIF
+
+        wo(:,:,1)=ozonecm(latitude_deg, paprs,read_climoz,rjour=zzz)
+      ELSE
+        ro3i = int((days_elapsed + jh_cur - jh_1jan) / year_len * 360.) + 1   
+        if (ro3i == 361) ro3i = 360
+        if (read_climoz == 1) then
+           call regr_pr_av(ncid_climoz, (/"tro3"/), julien=ro3i, &
+                press_in_edg=press_climoz, paprs=paprs, v3=wo)
+        else
+           ! read_climoz == 2
+           call regr_pr_av(ncid_climoz, (/"tro3         ", &
+                "tro3_daylight"/), julien=ro3i, press_in_edg=press_climoz, &
+                paprs=paprs, v3=wo)
+        end if
+        ! Convert from mole fraction of ozone to column density of ozone in a
+        ! cell, in kDU:
+        forall (l = 1: read_climoz) wo(:, :, l) = wo(:, :, l) * rmo3 / rmd &
+             * zmasse / dobson_u / 1e3
+        ! (By regridding ozone values for LMDZ only once every 360th of
+        ! year, we have already neglected the variation of pressure in one
+        ! 360th of year. So do not recompute "wo" at each time step even if
+        ! "zmasse" changes a little.)
+     
+      ENDIF
     ENDIF
     !