diff --git a/source/model/lake.f90 b/source/model/lake.f90
index 72b3311e0539b0486675e5095efe49a846b56f92..269301d717d3e37fd4c2ed7ee089708508d589b0 100644
--- a/source/model/lake.f90
+++ b/source/model/lake.f90
@@ -932,6 +932,9 @@ endif
! Calculation of the whole temperature profile
call T_SOLVER(ix,iy,dnx,dny,year,month,day,hour,phi, &
& RadWater, RadIce, SurfTemp, fetch, dt)
+if (soilswitch%par == 1) then
+ call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
+endif
! Diagnostic calculations
do i = 1, M
@@ -958,7 +961,6 @@ call OXYGEN_PRODCONS(gs, gsp, wst, bathymsoil, gas, area_int, area_half, ddz05,
& prodox, resp, bod, sod, sodbot)
if (soilswitch%par == 1) then
- call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
! Calculation of methane in all soil columns, except for the lowest one
if (multsoil) then
call SOILCOLSTEMP(gs,gsp,dt,ls,ftot,ch4_pres_atm,ddz,ddzi,zsoilcols, &
@@ -987,7 +989,7 @@ if (soilswitch%par == 1) then
gs%isoilcol = nsoilcols
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,nsoilcols),rosoil,fbbleflx_ch4_sum, &
- & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals2, &
+ & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals1, &
& rootss,rosdry,por,veg,0._ireals,TgrAnn, methgenmh, &
& ddz,ddz05,wst, lammeth, h1, ls, dhw, dhw0, .true., .false., &
& lsm, bathymwater, &
@@ -1316,6 +1318,9 @@ if (flag_snow == 1) then
! Calculation of the whole temperature profile
call T_SOLVER(ix,iy,dnx,dny,year,month,day,hour,phi, &
& RadWater, RadIce, SurfTemp, fetch, dt)
+ if (soilswitch%par == 1) then
+ call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
+ endif
!Assuming all dissolved species diffuse with the same molecular diffusivity as salinity
lamsal(:) = (lamw(:) - lamw0)/cw_m_row0 * alsal + lamsal0
@@ -1328,7 +1333,6 @@ if (flag_snow == 1) then
call S_DIFF(dt,ls,salice)
if (saltice%par == 1) call UPDATE_ICE_SALINITY(dt,ls,wst)
if (soilswitch%par == 1) then
- call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
! Calculation of methane in all soil columns, except for the lowest one
if (multsoil) then
call SOILCOLSTEMP(gs,gsp,dt,ls,ftot,ch4_pres_atm,ddz,ddzi,zsoilcols, &
@@ -1366,7 +1370,7 @@ if (flag_snow == 1) then
gs%isoilcol = nsoilcols
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,nsoilcols),rosoil,fbbleflx_ch4_sum, &
- & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals2, &
+ & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals1, &
& rootss,rosdry,por,veg,0._ireals,TgrAnn, methgenmh, &
& ddz,ddz05,wst, lammeth, h1, ls, dhw, dhw0, .true., .false. ,&
& lsm, bathymwater, &
@@ -1583,6 +1587,9 @@ else
! Calculation of the whole temperature profile
call T_SOLVER(ix,iy,dnx,dny,year,month,day,hour,phi, &
& RadWater, RadIce, SurfTemp, fetch, dt)
+ if (soilswitch%par == 1) then
+ call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
+ endif
!Assuming all dissolved species diffuse with the same molecular diffusivity as salinity
lamsal(:) = (lamw(:) - lamw0)/cw_m_row0 * alsal + lamsal0
@@ -1603,7 +1610,6 @@ else
& prodox, resp, bod, sod, sodbot)
if (soilswitch%par == 1) then
- call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
! Calculation of temperature in all soil columns, except for the lowest one
if (multsoil) then
call SOILCOLSTEMP(gs,gsp,dt,ls,ftot,ch4_pres_atm,ddz,ddzi,zsoilcols, &
@@ -1630,7 +1636,7 @@ else
gs%isoilcol = nsoilcols
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,nsoilcols),rosoil,fbbleflx_ch4_sum, &
- & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals2, &
+ & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals1, &
& rootss,rosdry,por,veg,0._ireals,TgrAnn, methgenmh, &
& ddz,ddz05,wst, lammeth, h1, ls, dhw, dhw0, .true., .false., &
& lsm, bathymwater, &
@@ -1813,13 +1819,13 @@ if (flag_snow == 1) then
! Temperature profile calculation
call T_SOLVER(ix,iy,dnx,dny,year,month,day,hour,phi, &
& RadWater, RadIce, SurfTemp, fetch, dt)
+ if (soilswitch%par == 1) then
+ call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
+ endif
! Salinity profile
call S_DIFF(dt,ls,salice)
if (saltice%par == 1) call UPDATE_ICE_SALINITY(dt,ls,wst)
- if (soilswitch%par == 1) then
- call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
- endif
call SNOWTEMP(ix,iy,dnx,dny,year,month,day,hour,snowmass, &
& snowmass_init,a,b,c,d,Temp,phi,fetch,dt)
if (soilswitch%par == 1) then
@@ -1827,7 +1833,7 @@ if (flag_snow == 1) then
gs%isoilcol = nsoilcols
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,nsoilcols),rosoil,fbbleflx_ch4_sum, &
- & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals2, &
+ & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals1, &
& rootss,rosdry,por,veg,0._ireals,TgrAnn, methgenmh, &
& ddz,ddz05,wst, lammeth, h1, ls, dhw, dhw0, .true., .false., &
& lsm, bathymwater, &
@@ -1898,17 +1904,19 @@ else
& gsp,soilflux(1,ix,iy),.false.,lsh)
call T_SOLVER(ix,iy,dnx,dny,year,month,day,hour,phi, &
& RadWater, RadIce, SurfTemp, fetch, dt)
+ if (soilswitch%par == 1) then
+ call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
+ endif
! Salinity profile
call S_DIFF(dt,ls,salice)
if (saltice%par == 1) call UPDATE_ICE_SALINITY(dt,ls,wst)
if (soilswitch%par == 1) then
- call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
call BUBBLE_BLOCK
gs%isoilcol = nsoilcols
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,nsoilcols),rosoil,fbbleflx_ch4_sum, &
- & fbbleflx_ch4(0,nsoilcols), wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals2, &
+ & fbbleflx_ch4(0,nsoilcols), wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals1, &
& rootss,rosdry,por,veg,0._ireals,TgrAnn, methgenmh, &
& ddz,ddz05,wst, lammeth, h1, ls, dhw, dhw0, .true., .false., &
& lsm, bathymwater, &
@@ -1954,12 +1962,12 @@ if4: IF (layer_case==4) THEN
endif
call T_SOLVER(ix,iy,dnx,dny,year,month,day,hour,phi, &
& RadWater, RadIce, SurfTemp, fetch, dt)
-
- ! Salinity profile
- call S_DIFF(dt,ls,salice)
if (soilswitch%par == 1) then
call SOILFORLAKE(dt,a,b,c,d,nsoilcols)
endif
+
+ ! Salinity profile
+ call S_DIFF(dt,ls,salice)
call SNOWTEMP(ix,iy,dnx,dny,year,month,day,hour,snowmass, &
& snowmass_init,a,b,c,d,Temp,phi,fetch,dt)
if (soilswitch%par == 1) then
@@ -1967,7 +1975,7 @@ if4: IF (layer_case==4) THEN
gs%isoilcol = nsoilcols
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,nsoilcols),rosoil,fbbleflx_ch4_sum, &
- & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals2, &
+ & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals1, &
& rootss,rosdry,por,veg,0._ireals,TgrAnn, methgenmh, &
& ddz,ddz05,wst, lammeth, h1, ls, dhw, dhw0, .true., .false., &
& lsm, bathymwater, &
@@ -2005,7 +2013,7 @@ if4: IF (layer_case==4) THEN
call BUBBLE_BLOCK
gs%isoilcol = nsoilcols
call METHANE(gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,nsoilcols),rosoil,fbbleflx_ch4_sum, &
- & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals2, &
+ & fbbleflx_ch4(0,nsoilcols),wl4(1,nsoilcols),wi2(1,nsoilcols),wa,Sals1, &
& rootss,rosdry,por,veg,0._ireals,TgrAnn, methgenmh, &
& ddz,ddz05,wst, lammeth, h1, ls, dhw, dhw0, .true., .false., &
& lsm, bathymwater, &
diff --git a/source/model/methane_mod.f90 b/source/model/methane_mod.f90
index b71925ae567380f036139e3d3bc1a75918c5fd07..c5ee1995f923b0ca0e80bbd75d947367d345d9b2 100644
--- a/source/model/methane_mod.f90
+++ b/source/model/methane_mod.f90
@@ -639,8 +639,6 @@ if (deepestsoil) then
call PROGONKA (a,b,c,f,y,1,gs%M+1)
y(1:gs%M+1) = max(y(1:gs%M+1),0.e0_ireals)
qwater(1:gs%M+1,2) = y(1:gs%M+1)
- !print*, lammeth(1), lammeth(gs%M)
- !print*, y(1), y(gs%M+1)
else
bottom_layer_if : if (ls%ls1 > 0) then
! Case water, deepice and soil; upper ice and snow are allowed
diff --git a/source/model/soil_mod.f90 b/source/model/soil_mod.f90
index b13e7069311296583f698881c11086c3808fcdfc..275f0831b8e08abac7966c81c013aa718351a703 100644
--- a/source/model/soil_mod.f90
+++ b/source/model/soil_mod.f90
@@ -478,7 +478,7 @@ enddo
! STEP 3 OF SPLITTING-UP METHOD : PHASE PROCESSES
mhsep = (1. + tricemethhydr%par*alphamh)
20 c2: do i = 1, ns
- work = MELTINGPOINT(Sals2(i,is)/wl3(i),pressoil(i),tricemethhydr%par)
+ work = MELTINGPOINT(Sals1(i,is)/wl3(i),pressoil(i),tricemethhydr%par)
if (wi1(i,is) > 0 .and. Tsoil2(i) > work + 0.01) then
Potphenergy = (Tsoil2(i) - (work + 0.01)) * &
& csoil(i)*rosoil(i)*dzss(i)
@@ -534,7 +534,7 @@ enddo c2
max1 = 0.
do i = 1, ns
if (Tsoil3(i,is) < &
- & MELTINGPOINT(Sals2(i,is)/wl4(i,is),pressoil(i),tricemethhydr%par) - 0.01 &
+ & MELTINGPOINT(Sals1(i,is)/wl4(i,is),pressoil(i),tricemethhydr%par) - 0.01 &
&.and. abs(wl4(i,is) - UNFRWAT(Tsoil3(i,is),i)) > max1) then
max1 = abs(wl4(i,is) - UNFRWAT(Tsoil3(i,is),i))
endif
@@ -807,7 +807,7 @@ use PHYS_FUNC, only : &
use ARRAYS_BATHYM, only : bathym, dhw, dhw0, layers_type
use ARRAYS_SOIL, only : csoil, rosoil, lamsoil, Tsoil3, Tsoil2, Tsoil1, &
-& wl4,wi2,wa,Sals2,rosdry,por,lsm,dzs,zsoil,rosoil
+& wl4,wi2,wa,Sals1,rosdry,por,lsm,dzs,zsoil,rosoil
use ARRAYS_GRID, only : ddz05, gridsize_type, gridspacing_type
use ARRAYS_METHANE, only : veg,TgrAnn,methgenmh,qwater,lammeth, &
& fplant,febul0, fdiff_lake_surf, &
@@ -942,7 +942,7 @@ if (contr(2) == 1) then
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,is),rosoil, &
& work, work, &
- & wl4(1,is),wi2(1,is),wa,Sals2, &
+ & wl4(1,is),wi2(1,is),wa,Sals1, &
& rootss,rosdry,por,veg,qsoilsurf(is)*por(1),TgrAnn, methgenmh, &
& ddz,ddz05, wst, lammeth, &
& gsp%z_half(bathymsoil(is)%icent), & !0.5*(zsoilcols(is) + zsoilcols(is+1)), &
@@ -979,7 +979,7 @@ if (contr(2) == 1) then
call METHANE &
& (gas,pressure,wind10,ch4_pres_atm,zsoil,Tsoil3(1,is),rosoil, &
& work, work, &
- & wl4(1,is),wi2(1,is),wa,Sals2, &
+ & wl4(1,is),wi2(1,is),wa,Sals1, &
& rootss,rosdry,por,veg,fdiffbot(is),TgrAnn, methgenmh, &
& ddz,ddz05, wst, lammeth, &
& gsp%z_half(bathymsoil(is)%icent), & !0.5*(zsoilcols(is) + zsoilcols(is+1)), &