Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Awicm3 ice thermo cpl fixes #93

Merged
merged 6 commits into from
Apr 12, 2021
Merged

Awicm3 ice thermo cpl fixes #93

Changes from all commits
Commits
Show all changes
6 commits
Select commit Hold shift + click to select a range
8ef69a2
including hotfix for sea ice initial temperature in previously open w…
JanStreffing Mar 3, 2021
d3924fb
test setup for RHO1
JanStreffing Mar 3, 2021
d248fae
test for FRE1
JanStreffing Mar 3, 2021
8f50ef1
test FRE2
JanStreffing Mar 3, 2021
9ddc871
adding SH/NH lead_closing parameter
JanStreffing Mar 3, 2021
e3c08ef
removing faulty density scaling
JanStreffing Apr 10, 2021
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
42 changes: 23 additions & 19 deletions src/ice_thermo_cpl.F90
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -132,17 +132,17 @@ subroutine thermodynamics(mesh)
rsf = 0._WP
end if

!#if defined (__oifs)
! !---- different lead closing parameter for NH and SH
! call r2g(geolon, geolat, coord_nod2d(1,inod), coord_nod2d(2,inod))
! if (geolat.lt.0.) then
! h0min = 0.75
! h0max = 1.0
! else
! h0min = 0.5
! h0max = 0.75
! endif
!#endif /* (__oifs) */
#if defined (__oifs)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

We had before the situation, when people where trying to optimize for h0, available in the namelist, but then hardcoded downstream. Will we get into similar trap here? Should we maybe expose those parameters explicitly through namelist.ice and comment that they only used in a coupled setting?

Copy link
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

great job, great bug fixes! technically it can be merged as is. @koldunovn shall we introduce namelist.cpl or a section &cpl in one of the namelists for FESOM?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Maybe not for this PR, but later I think it's not a bad idea. If you OK with this PR I would merge as is.

Copy link
Collaborator Author

@JanStreffing JanStreffing Apr 12, 2021 •
edited
Loading

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I would be ok with that. And to have the namelist.cpl be a separate development.

!---- different lead closing parameter for NH and SH
call r2g(geolon, geolat, coord_nod2d(1,inod), coord_nod2d(2,inod))
if (geolat.lt.0.) then
h0min = 1.0
h0max = 1.0
else
h0min = 0.3
h0max = 0.3
endif
#endif /* (__oifs) */

call ice_growth
#if defined (__oifs)
Expand All @@ -154,7 +154,8 @@ subroutine thermodynamics(mesh)
call ice_surftemp(max(h/(max(A,Aimin)),0.05),hsn/(max(A,Aimin)),a2ihf,t)
ice_temp(inod) = t
else
ice_temp(inod) = 275.15_WP
! Freezing temp of saltwater in K
ice_temp(inod) = -0.0575_WP*S_oc_array(inod) + 1.7105e-3_WP*sqrt(S_oc_array(inod)**3) -2.155e-4_WP*(S_oc_array(inod)**2)+273.15_WP
endif
call ice_albedo(h,hsn,t,alb)
ice_alb(inod) = alb
Expand Down Expand Up @@ -473,18 +474,21 @@ subroutine ice_surftemp(h,hsn,a2ihf,t)
real(kind=WP) zcprosn
!---- local parameters
real(kind=WP), parameter :: dice = 0.05_WP ! ECHAM6's thickness for top ice "layer"
real(kind=WP), parameter :: ctfreez = 271.38_WP ! ECHAM6's temperature at which sea starts freezing/melting
!---- freezing temperature of sea-water [K]
real(kind=WP) :: TFrezs

!---- compute freezing temperature of sea-water from salinity
TFrezs = -0.0575_WP*S_oc + 1.7105e-3_WP*sqrt(S_oc**3) - 2.155e-4_WP*(S_oc**2)+273.15

snicecond = con/consn*rhowat/rhosno ! equivalence fraction thickness of ice/snow
zsniced=h+snicecond*hsn ! Ice + Snow-Ice-equivalent thickness [m]
zicefl=con*ctfreez/zsniced ! Conductive heat flux through sea ice [W/m²]
snicecond = con/consn ! equivalence fraction thickness of ice/snow
zsniced=h+snicecond*hsn ! Ice + Snow-Ice-equivalent thickness [m]
zicefl=con*TFrezs/zsniced ! Conductive heat flux through sea ice [W/m²]
hcapice=rhoice*cpice*dice ! heat capacity of upper 0.05 cm sea ice layer [J/(m²K)]
zcpdt=hcapice/dt ! Energy required to change temperature of top ice "layer" [J/(sm²K)]
zcprosn=rhowat*cpsno/dt ! Specific Energy required to change temperature of 1m snow on ice [J/(sm³K)]
zcprosn=rhosno*cpsno/dt ! Specific Energy required to change temperature of 1m snow on ice [J/(sm³K)]
zcpdte=zcpdt+zcprosn*hsn ! Combined Energy required to change temperature of snow + 0.05m of upper ice
t=(zcpdte*t+a2ihf+zicefl)/(zcpdte+con/zsniced) ! New sea ice surf temp [K]
t=min(ctfreez,t) ! Not warmer than freezing please!
t=min(TFrezs,t) ! Not warmer than freezing please!
end subroutine ice_surftemp

subroutine ice_albedo(h,hsn,t,alb)
Expand All @@ -505,7 +509,7 @@ subroutine ice_albedo(h,hsn,t,alb)
! set albedo
! ice and snow, freezing and melting conditions are distinguished
if (h>0.0_WP) then
if (t<273.14_WP) then ! freezing condition
if (t<273.15_WP) then ! freezing condition
if (hsn.gt.0.0_WP) then ! snow cover present
alb=albsn
else ! no snow cover
Expand Down