Changeset 2430

Timestamp:

Nov 6, 2020, 2:26:13 PM (4 years ago)

Author:

emillour

Message:

Venus GCM:
Improve hgardfou by adding checks for NaN in the tests.
EM

File:

• trunk/LMDZ.VENUS/libf/phyvenus/hgardfou.F (modified) (3 diffs)

trunk/LMDZ.VENUS/libf/phyvenus/hgardfou.F

@@ -4 +4 @@
       SUBROUTINE hgardfou (t,tsol,text)
 
-      use dimphy
+      use dimphy, only: klon,klev
       IMPLICIT none
 c======================================================================
-c Verifier la temperature
+c Check temperature
 c======================================================================
-#include "YOMCST.h"
-      REAL t(klon,klev), tsol(klon)
-      CHARACTER*(*) text
+!#include "YOMCST.h"
+      REAL,INTENT(IN) :: t(klon,klev), tsol(klon)
+      CHARACTER(len=*),INTENT(in):: text
 C
-      INTEGER i, k
-      REAL zt(klon)
-      INTEGER jadrs(klon), jbad
-      LOGICAL ok
+      INTEGER :: i, k
+      REAL :: zt(klon) ! to locally store temperature fields
+      REAL,PARAMETER :: tmin=15.0   ! Minimum temperature
+      REAL,PARAMETER :: tmax=1200.0 ! Maximum temperature
+      INTEGER :: jbad ! number of problematic points
+      INTEGER :: jadrs(klon) ! stored index of problematic points
+      LOGICAL :: ok ! .true. as long as everything OK
+      CHARACTER(len=100) :: textout
 c
-      LOGICAL firstcall
-      SAVE firstcall
-      DATA firstcall /.TRUE./
+      LOGICAL,SAVE :: firstcall=.TRUE.
+
       IF (firstcall) THEN
-         PRINT*, 'hgardfou garantit la temperature dans [15,1200] K'
+         PRINT*, 'hgardfou checks if temperature is in [15,1200] K'
          firstcall = .FALSE.
       ENDIF
 c
       ok = .TRUE.
+
+      !1. Atmospheric temperatures 
       DO k = 1, klev
          DO i = 1, klon
-!!!! MODIF GG to avoid crash after 78--> 95 extension!!
-!!!  WARNING: it has to be review/removed when the extension to the
-!!!  thermosphere is completed (physical processes and ionosphere added)
-          IF (k.LT.85) THEN
             zt(i) = t(i,k)
-          ENDIF
          ENDDO
+
+         ! Look for temperatures that are not numbers (NaN, Infinity, etc.)
+         jbad=0
+         DO i=1,klon
+           IF (zt(i).NE.zt(i)) THEN
+             jbad = jbad + 1
+             jadrs(jbad) = i
+           ENDIF
+         ENDDO
+         
+         IF (jbad .GT. 0) THEN
+           ok = .FALSE.
+           DO i = 1, jbad
+             PRINT *,'i,k,temperature =',jadrs(i),k,zt(jadrs(i))
+           ENDDO
+         ENDIF
+         
+         
+         ! Look for temperatures greater than tmax
 #ifdef CRAY
-         CALL WHENFGT(klon, zt, 1, 1200.0, jadrs, jbad)
+         CALL WHENFGT(klon, zt, 1, tmax, jadrs, jbad)
 #else
          jbad = 0
          DO i = 1, klon
-         IF (zt(i).GT.1200.0) THEN
+           IF (zt(i).GT.tmax) THEN
             jbad = jbad + 1
             jadrs(jbad) = i
-         ENDIF
+           ENDIF
          ENDDO
 #endif
@@ -53 +72 @@
            ENDDO
          ENDIF
+
+         ! Look for temperatures lower than tmin
 #ifdef CRAY
-         CALL WHENFLT(klon, zt, 1, 15.0, jadrs, jbad)
+         CALL WHENFLT(klon, zt, 1, tmin, jadrs, jbad)
 #else
          jbad = 0
          DO i = 1, klon
-         IF (zt(i).LT.15.0) THEN
+           IF (zt(i).LT.tmin) THEN
             jbad = jbad + 1
             jadrs(jbad) = i
-         ENDIF
+           ENDIF
          ENDDO
 #endif
@@ -70 +91 @@
            ENDDO
          ENDIF
+      ENDDO ! of DO k = 1, klev
+
+      !2. surface temperatures
+      DO i = 1, klon
+         zt(i) = tsol(i)
       ENDDO
-c
-         DO i = 1, klon
-            zt(i) = tsol(i)
-         ENDDO
+
+      ! Look for temperatures that are not numbers (NaN, Infinity, etc.)
+      jbad=0
+      DO i=1,klon
+        IF (zt(i).NE.zt(i)) THEN
+          jbad = jbad + 1
+          jadrs(jbad) = i
+        ENDIF
+      ENDDO
+         
+      IF (jbad .GT. 0) THEN
+        ok = .FALSE.
+        DO i = 1, jbad
+          PRINT *,'i,temperature =',jadrs(i),zt(jadrs(i))
+        ENDDO
+      ENDIF
+
+      ! Look for temperatures greater then tmax
 #ifdef CRAY
-         CALL WHENFGT(klon, zt, 1, 1200.0, jadrs, jbad)
+      CALL WHENFGT(klon, zt, 1, tmax, jadrs, jbad)
 #else
-         jbad = 0
-         DO i = 1, klon
-         IF (zt(i).GT.1200.0) THEN
-            jbad = jbad + 1
-            jadrs(jbad) = i
-         ENDIF
-         ENDDO
+      jbad = 0
+      DO i = 1, klon
+        IF (zt(i).GT.tmax) THEN
+           jbad = jbad + 1
+           jadrs(jbad) = i
+        ENDIF
+      ENDDO
 #endif
-         IF (jbad .GT. 0) THEN
-           ok = .FALSE.
-           DO i = 1, jbad
-             PRINT *,'i,temperature =',jadrs(i),zt(jadrs(i))
-           ENDDO
-         ENDIF
+      IF (jbad .GT. 0) THEN
+        ok = .FALSE.
+        DO i = 1, jbad
+          PRINT *,'i,temperature =',jadrs(i),zt(jadrs(i))
+        ENDDO
+      ENDIF
+
+      ! Look for temperatures lower than tmin
 #ifdef CRAY
-         CALL WHENFLT(klon, zt, 1, 20.0, jadrs, jbad)
+      CALL WHENFLT(klon, zt, 1, tmin, jadrs, jbad)
 #else
-         jbad = 0
-         DO i = 1, klon
-         IF (zt(i).LT.20.0) THEN
-            jbad = jbad + 1
-            jadrs(jbad) = i
-         ENDIF
-         ENDDO
+      jbad = 0
+      DO i = 1, klon
+        IF (zt(i).LT.tmin) THEN
+           jbad = jbad + 1
+           jadrs(jbad) = i
+        ENDIF
+      ENDDO
 #endif
-         IF (jbad .GT. 0) THEN
-           ok = .FALSE.
-           DO i = 1, jbad
-             PRINT *,'i,temperature =',jadrs(i),zt(jadrs(i))
-           ENDDO
-         ENDIF
+      IF (jbad .GT. 0) THEN
+        ok = .FALSE.
+        DO i = 1, jbad
+          PRINT *,'i,temperature =',jadrs(i),zt(jadrs(i))
+        ENDDO
+      ENDIF
 c
       IF (.NOT. ok) THEN
-         text='hgardfou s arrete '//text
-         CALL abort_physic("hgardfou", text, 1)
+         textout='hgardfou stops '//text
+         CALL abort_physic("hgardfou", textout, 1)
       ENDIF
 
-      RETURN
       END