Changeset 18adf60 in flexpart.git for src/advance.f90

Timestamp:

May 9, 2016, 9:35:01 AM (8 years ago)

Author:

Espen Sollum ATMOS <eso@…>

Branches:

master, 10.4.1_pesei, GFS_025, bugfixes+enhancements, dev, release-10, release-10.4.1, scaling-bug, univie

Children:

3b80e98, 93786a1

Parents:

62e65c7

Message:

Bug found in the MQUASILAG=1 option, xmass array accessed with out-of-bounds index. Quick fix implemented (see changes in timemanager and advance).

File:

• src/advance.f90 (modified) (3 diffs)

src/advance.f90

@@ -493 +493 @@
 
       if (mdomainfill.eq.0) then
-        do nsp=1,nspec
-          if (xmass(nrelpoint,nsp).gt.eps2) goto 887
-        end do
-887     nsp=min(nsp,nspec)
-!!$        if (density(nsp).gt.0.) &
-!!$             call get_settling(itime,xts,yts,zt,nsp,settling)    !old
-        if (density(nsp).gt.0.) &
-             call get_settling(itime,real(xt),real(yt),zt,nsp,settling)  !bugfix
+! ESO 05.2015  Changed this to fix MQUASILAG option, where nrelpoint is
+!              particle number and thus xmass array goes out of bounds
+!        do nsp=1,nspec
+!           if (xmass(nrelpoint,nsp).gt.eps2) goto 887
+!         end do
+! 887     nsp=min(nsp,nspec)
+        if (nspec.eq.1.and.density(1).gt.0.) then
+          call get_settling(itime,real(xt),real(yt),zt,nsp,settling)  !bugfix
+        end if
         w=w+settling
       endif
@@ -654 +655 @@
 
     if (mdomainfill.eq.0) then
-      do nsp=1,nspec
-        if (xmass(nrelpoint,nsp).gt.eps2) goto 888
-      end do
-888   nsp=min(nsp,nspec)
-!!$      if (density(nsp).gt.0.) &
-!!$           call get_settling(itime,xts,yts,zt,nsp,settling)    !old
-      if (density(nsp).gt.0.) &
-           call get_settling(itime,real(xt),real(yt),zt,nsp,settling)  !bugfix
+! ESO 05.2015  Changed this to fix MQUASILAG option, where nrelpoint is
+!              particle number and thus xmass array goes out of bounds
+
+!      do nsp=1,nspec
+!         if (xmass(nrelpoint,nsp).gt.eps2) goto 888
+!       end do
+! 888   nsp=min(nsp,nspec)
+!        if (density(nsp).gt.0.) then
+      if (nspec.eq.1.and.density(1).gt.0.) then
+        call get_settling(itime,real(xt),real(yt),zt,nsp,settling)  !bugfix
+      end if
       w=w+settling
     endif
@@ -858 +862 @@
 
   if (mdomainfill.eq.0) then
-    do nsp=1,nspec
-      if (xmass(nrelpoint,nsp).gt.eps2) goto 889
-    end do
-889   nsp=min(nsp,nspec)
-!!$    if (density(nsp).gt.0.) &
-!!$         call get_settling(itime+ldt,xts,yts,zt,nsp,settling)    !old
-    if (density(nsp).gt.0.) &
-         call get_settling(itime+ldt,real(xt),real(yt),zt,nsp,settling)  !bugfix
+! ESO 05.2015  Changed this to fix MQUASILAG option, where nrelpoint is
+!              particle number and thus xmass array goes out of bounds
+!    do nsp=1,nspec
+!       if (xmass(nrelpoint,nsp).gt.eps2) goto 889
+!     end do
+! 889   nsp=min(nsp,nspec)
+!      if (density(nsp).gt.0.) then
+    if (nspec.eq.1.and.density(1).gt.0.) then
+      call get_settling(itime+ldt,real(xt),real(yt),zt,nsp,settling)  !bugfix
+    end if
     w=w+settling
   endif