fix array-positions when useing merge for bombs (output vs running po…

…sition)
metno · Jan 26, 2025 · eab0c2f · eab0c2f
1 parent e861dd1
commit eab0c2f
Show file tree

Hide file tree

Showing 3 changed files with 20 additions and 17 deletions.
diff --git a/src/common/decay.f90 b/src/common/decay.f90
@@ -50,13 +50,14 @@ end subroutine decay
 subroutine decayDeps(tstep)
   USE snapfldML, only: depdry, depwet, accdry, accwet, &
     total_activity_released, total_activity_lost_domain, total_activity_lost_other
-  USE snapparML, only: ncomp, run_comp, def_comp
+  USE snapparML, only: ncomp, nocomp, run_comp, def_comp, output_component
   USE iso_fortran_env, only: real64
   USE releaseML, only: tpos_bomb
 
   real, intent(in) :: tstep
 
   integer :: m, mm
+  integer, allocatable :: mmo(:)
   real :: bomb_decay_rate, current_state, next_state
   logical, save :: prepare = .TRUE.
   logical, save :: has_bomb_decay = .FALSE.
@@ -103,17 +104,19 @@ subroutine decayDeps(tstep)
     total_steps = total_steps + tstep
   end if
 
-  do m=1,ncomp
-    mm = run_comp(m)%to_defined
-    if(def_comp(mm)%kdecay >= 1) then
-      depdry(:,:,m) = depdry(:,:,m)*def_comp(mm)%decayrate
-      depwet(:,:,m) = depwet(:,:,m)*def_comp(mm)%decayrate
-      accdry(:,:,m) = accdry(:,:,m)*def_comp(mm)%decayrate
-      accwet(:,:,m) = accwet(:,:,m)*def_comp(mm)%decayrate
-
-      total_activity_released(m) = total_activity_released(m)*def_comp(mm)%decayrate
-      total_activity_lost_domain(m) = total_activity_lost_domain(m)*def_comp(mm)%decayrate
-      total_activity_lost_other(m) = total_activity_lost_other(m)*def_comp(mm)%decayrate
+  do m=1,nocomp
+    mmo = output_component(m)%to_defined
+    ! decay-rates of merged components must be similar, otherwise, deday of output-fields
+    ! does not work well enough
+    if(def_comp(mmo(1))%kdecay >= 1) then
+      depdry(:,:,m) = depdry(:,:,m)*def_comp(mmo(1))%decayrate
+      depwet(:,:,m) = depwet(:,:,m)*def_comp(mmo(1))%decayrate
+      accdry(:,:,m) = accdry(:,:,m)*def_comp(mmo(1))%decayrate
+      accwet(:,:,m) = accwet(:,:,m)*def_comp(mmo(1))%decayrate
+
+      total_activity_released(m) = total_activity_released(m)*def_comp(mmo(1))%decayrate
+      total_activity_lost_domain(m) = total_activity_lost_domain(m)*def_comp(mmo(1))%decayrate
+      total_activity_lost_other(m) = total_activity_lost_other(m)*def_comp(mmo(1))%decayrate
     endif
   enddo
   return

diff --git a/src/common/drydep.f90 b/src/common/drydep.f90
@@ -48,7 +48,7 @@ subroutine drydep1(part)
       dep = part%scale_rad(1.0 - def_comp(m)%drydeprat)
       i = hres_pos(part%x)
       j = hres_pos(part%y)
-      mm = def_comp(m)%to_running
+      mm = def_comp(m)%to_output
     !$OMP atomic
       depdry(i,j,mm) = depdry(i,j,mm) + dble(dep)
     end if
@@ -98,7 +98,7 @@ subroutine drydep2(tstep, part)
     dep = part%scale_rad(1.0 - deprate)
     i = hres_pos(part%x)
     j = hres_pos(part%y)
-    mm = def_comp(m)%to_running
+    mm = def_comp(m)%to_output
   !$OMP atomic
     depdry(i,j,mm) = depdry(i,j,mm) + dble(dep)
   end if

diff --git a/src/common/wetdep.f90 b/src/common/wetdep.f90
@@ -44,7 +44,7 @@ subroutine wetdep2(depwet, tstep, part, pextra)
 !> uses the precipitation at the particle position
     type(extraParticle), intent(in) :: pextra
 
-    integer :: m, i, j, mm
+    integer :: m, i, j, mm, mo
     real :: precint, deprate, dep, q
 
     m = part%icomp
@@ -63,10 +63,10 @@ subroutine wetdep2(depwet, tstep, part, pextra)
 
       i = hres_pos(part%x)
       j = hres_pos(part%y)
-      mm = def_comp(m)%to_running
+      mo = def_comp(m)%to_output
 
       !$OMP atomic
-      depwet(i, j, mm) = depwet(i, j, mm) + dble(dep)
+      depwet(i, j, mo) = depwet(i, j, mo) + dble(dep)
     end if
   end subroutine wetdep2