Changeset 4332 for palm/trunk/SOURCE/synthetic_turbulence_generator_mod.f90

Timestamp:

Dec 10, 2019 7:44:12 PM (5 years ago)

Author:

suehring

Message:

Synthetic turbulence: Limit initial velocity seeds in restart runs, if not the seed calculation may become unstable.

File:

• palm/trunk/SOURCE/synthetic_turbulence_generator_mod.f90 (modified) (3 diffs)

palm/trunk/SOURCE/synthetic_turbulence_generator_mod.f90

@@ -25 +25 @@
 ! -----------------
 ! $Id$
+! Limit initial velocity seeds in restart runs, if not the seed calculation 
+! may become unstable. Further, minor bugfix in initial velocity seed 
+! calculation.
+! 
+! 4329 2019-12-10 15:46:36Z motisi
 ! Renamed wall_flags_0 to wall_flags_static_0
 ! 
@@ -903 +908 @@
           DO  j = nysg, nyng
              DO  k = nzb, nzt+1
-                IF  ( a11(k) /= 0.0_wp )  THEN
+                IF  ( a11(k) > 10E-8_wp )  THEN
                    fu_yz(k,j) = ( u(k,j,i) - u_init(k) ) / a11(k)
                 ELSE
-                   fu_yz(k,j) = 0.0_wp
+                   fu_yz(k,j) = 10E-8_wp
                 ENDIF
 
-                IF  ( a22(k) /= 0.0_wp )  THEN
+                IF  ( a22(k) > 10E-8_wp )  THEN
                    fv_yz(k,j) = ( v(k,j,i) -                                  &
                                   a21(k) * fu_yz(k,j) - v_init(k) ) / a22(k)
                 ELSE
-                   fv_yz(k,j) = 0.0_wp
+                   fv_yz(k,j) = 10E-8_wp
                 ENDIF
 
-                IF  ( a33(k) /= 0.0_wp )  THEN
+                IF  ( a33(k) > 10E-8_wp )  THEN
                    fw_yz(k,j) = ( w(k,j,i) -                                   &
                                   a31(k) * fu_yz(k,j) - a32(k) *               &
                                   fv_yz(k,j) ) / a33(k)
                 ELSE
-                   fw_yz = 0._wp
+                   fw_yz(k,j) = 10E-8_wp
                 ENDIF
-
              ENDDO
           ENDDO
@@ -936 +940 @@
              DO  k = nzb, nzt+1
 
-                IF  ( a11(k) /= 0.0_wp )  THEN
+                IF  ( a11(k) > 10E-8_wp )  THEN
                    fu_xz(k,i) = ( u(k,j,i) - u_init(k) ) / a11(k)
                 ELSE
-                   fu_xz(k,i) = 0.0_wp
+                   fu_xz(k,i) = 10E-8_wp
                 ENDIF
 
-                IF  ( a22(k) /= 0.0_wp )  THEN
+                IF  ( a22(k) > 10E-8_wp )  THEN
                    fv_xz(k,i) = ( v(k,j,i) -                                   &
                                   a21(k) * fu_xz(k,i) - v_init(k) ) / a22(k)
                 ELSE
-                   fv_xz(k,i) = 0.0_wp
+                   fv_xz(k,i) = 10E-8_wp
                 ENDIF
 
-                IF  ( a33(k) /= 0.0_wp )  THEN
+                IF  ( a33(k) > 10E-8_wp )  THEN
                    fw_xz(k,i) = ( w(k,j,i) -                                   &
                                   a31(k) * fu_xz(k,i) -                        &
                                   a32(k) * fv_xz(k,i) ) / a33(k)
                 ELSE
-                   fw_xz = 0._wp
+                   fw_xz(k,i) = 10E-8_wp
                 ENDIF