Changeset 3528

Timestamp:

Nov 15, 2018 7:07:11 PM (6 years ago)

Author:

suehring

Message:

Bugfix in raytracing - add an epsilon value to overcome precision-related errors.

File:

• palm/trunk/SOURCE/radiation_model_mod.f90 (modified) (6 diffs)

palm/trunk/SOURCE/radiation_model_mod.f90

@@ -28 +28 @@
 ! -----------------
 ! $Id$
+! Add an epsilon value to compare values in if statement to fix possible 
+! precsion related errors in raytrace routines.
+! 
+! 3524 2018-11-14 13:36:44Z raasch
 ! missing cpp-directives added
 ! 
@@ -7084 +7088 @@
         INTEGER(iwp)                           :: ip           !< number of processor where gridbox reside
         INTEGER(iwp)                           :: ig           !< 1D index of gridbox in global 2D array
-        REAL(wp)                               :: lad_s_target !< recieved lad_s of particular grid box
+        
+        REAL(wp)                               :: eps = 1E-10_wp !< epsilon for value comparison
+        REAL(wp)                               :: lad_s_target   !< recieved lad_s of particular grid box
 
 !
@@ -7167 +7173 @@
             ENDIF
 
-            IF ( nextdist >= distance ) EXIT
+            IF ( ABS(distance - nextdist) < eps )  EXIT
             lastdist = nextdist
             dimnext(seldim) = dimnext(seldim) + dimdelta(seldim)
@@ -7286 +7292 @@
       INTEGER(iwp)                           ::  nly          !< maximum  plant canopy height 
       INTEGER(iwp)                           ::  ntrack
-      REAL(wp)                               ::  zbottom, ztop !< urban surface boundary in real numbers
-      REAL(wp)                               ::  zorig        !< z coordinate of ray column entry
-      REAL(wp)                               ::  zexit        !< z coordinate of ray column exit
-      REAL(wp)                               ::  qdist        !< ratio of real distance to z coord difference
-      REAL(wp)                               ::  dxxyy        !< square of real horizontal distance
-      REAL(wp)                               ::  curtrans     !< transparency of current PC box crossing
+      
       INTEGER(iwp)                           ::  zb0
       INTEGER(iwp)                           ::  zb1
@@ -7304 +7305 @@
       INTEGER(MPI_ADDRESS_KIND)              ::  wdisp        !< RMA window displacement
 #endif
+      
+      REAL(wp)                               ::  eps = 1E-10_wp !< epsilon for value comparison
+      REAL(wp)                               ::  zbottom, ztop !< urban surface boundary in real numbers
+      REAL(wp)                               ::  zorig         !< z coordinate of ray column entry
+      REAL(wp)                               ::  zexit         !< z coordinate of ray column exit
+      REAL(wp)                               ::  qdist         !< ratio of real distance to z coord difference
+      REAL(wp)                               ::  dxxyy         !< square of real horizontal distance
+      REAL(wp)                               ::  curtrans      !< transparency of current PC box crossing
+      
+
       
       yxorigin(:) = origin(2:3)
@@ -7422 +7433 @@
          ENDIF
 
-         IF ( nextdist >= distance )  EXIT
+         IF ( ABS(distance - nextdist) < eps )  EXIT
 
          IF ( rad_angular_discretization  .AND.  calc_svf )  THEN