Changeset 5184a7c in flexpart.git for src/advance.f90

Timestamp:

Jun 20, 2017, 9:09:35 AM (7 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:

b5127f9

Parents:

a76d954

Message:

Enable settling with multiple species if from separate releases

File:

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

src/advance.f90

@@ -535 +535 @@
 
       if (mdomainfill.eq.0) then
-! 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,nspec,settling)  !bugfix
+        if (lsettling) then
+          do nsp=1,nspec
+            if (xmass(nrelpoint,nsp).gt.eps2) exit
+          end do
+          if (nsp.gt.nspec) then
+  ! This should never happen          
+            write(*,*) 'advance.f90: ERROR: could not find releasepoint'
+            stop
+          end if
+          if (density(nsp).gt.0.) then
+            call get_settling(itime,real(xt),real(yt),zt,nsp,settling)  !bugfix
+            w=w+settling
+          end if
         end if
-        w=w+settling
       endif
 
@@ -700 +704 @@
   !*************************************************************************
 
-
-
-    if (mdomainfill.eq.0) then
-! 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,nspec,settling)  !bugfix
+  if (mdomainfill.eq.0) then
+    if (lsettling) then
+      do nsp=1,nspec
+        if (xmass(nrelpoint,nsp).gt.eps2) exit
+      end do
+      if (nsp.gt.nspec) then
+  ! This should never happen          
+        write(*,*) 'advance.f90: ERROR: could not find releasepoint'
+        stop
       end if
-      w=w+settling
+      if (density(nsp).gt.0.) then
+        call get_settling(itime,real(xt),real(yt),zt,nsp,settling)  !bugfix
+        w=w+settling
+      end if
     endif
+  end if
+
 
   ! Calculate position at time step itime+lsynctime
@@ -910 +915 @@
 
   if (mdomainfill.eq.0) then
-! 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,nspec,settling)  !bugfix
-    end if
-    w=w+settling
-  endif
+    if (lsettling) then
+      do nsp=1,nspec
+        if (xmass(nrelpoint,nsp).gt.eps2) exit
+      end do
+      if (nsp.gt.nspec) then
+  ! This should never happen          
+        write(*,*) 'advance.f90: ERROR: could not find releasepoint'
+        stop
+      end if
+      if (density(nsp).gt.0.) then
+        call get_settling(itime+ldt,real(xt),real(yt),zt,nsp,settling) !bugfix
+        w=w+settling
+      end if
+    endif
+  end if