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.

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

Implement sea-ice freeboard cost #821

Open
wants to merge 1 commit into
base: master
Choose a base branch
from
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Jump to
Jump to file
Failed to load files.
Diff view
Diff view
3 changes: 3 additions & 0 deletions pkg/ecco/cost_gencost_customize.F
Original file line number Diff line number Diff line change
Expand Up @@ -221,6 +221,9 @@ subroutine cost_gencost_customize( mythid )
elseif (gencost_barfile(k)(1:9).EQ.'m_sihsnow') then
gencost_modfld(i,j,bi,bj,k) =
& hsnow(i,j,bi,bj)*maskC(i,j,1,bi,bj)
elseif (gencost_barfile(k)(1:11).EQ.'m_freeboard') then
gencost_modfld(i,j,bi,bj,k) =
& gencost_storefld(i,j,bi,bj,k)*maskC(i,j,1,bi,bj)
#endif
#ifdef ALLOW_GENCOST3D
elseif (gencost_barfile(k)(1:7).EQ.'m_theta') then
Expand Down
41 changes: 39 additions & 2 deletions pkg/ecco/ecco_phys.F
Original file line number Diff line number Diff line change
@@ -1,4 +1,7 @@
#include "ECCO_OPTIONS.h"
#ifdef ALLOW_SEAICE
# include "SEAICE_OPTIONS.h"
#endif
#ifdef ALLOW_SHELFICE
# include "SHELFICE_OPTIONS.h"
#endif
Expand Down Expand Up @@ -30,6 +33,11 @@ SUBROUTINE ECCO_PHYS( myTime, myIter, myThid )
# include "ECCO_SIZE.h"
# include "ECCO.h"
#endif
#ifdef ALLOW_SEAICE
# include "SEAICE_SIZE.h"
# include "SEAICE.h"
# include "SEAICE_PARAMS.h"
#endif
#ifdef ALLOW_PTRACERS
# include "PTRACERS_SIZE.h"
# include "PTRACERS_FIELDS.h"
Expand Down Expand Up @@ -82,16 +90,23 @@ SUBROUTINE ECCO_PHYS( myTime, myIter, myThid )
#ifdef ALLOW_PSBAR_STERIC
_RL VOLsumTile(nSx,nSy), RHOsumTile(nSx,nSy)
_RL VOLsumGlob_1, RHOsumGlob_1
c CHARACTER*(MAX_LEN_MBUF) msgBuf
C CHARACTER*(MAX_LEN_MBUF) msgBuf
#endif
C Mload :: total mass load (kg/m**2)
c _RL Mload(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
C _RL Mload(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
#ifdef ALLOW_SEAICE
_RL tmprecip_area, tmphsnowact, tmpheffact
_RL area_reg_sq
#endif
#if ( defined ALLOW_AUTODIFF_TAMC && defined ALLOW_GENCOST_CONTRIBUTION )
C ikey :: tape key (depends on kgen)
INTEGER ikey
#endif
CEOP

#ifdef ALLOW_SEAICE
area_reg_sq = SEAICE_area_reg * SEAICE_area_reg
#endif
tmpfac = recip_rhoConst*recip_gravity
sIceLoadFacLoc = zeroRL
IF ( useRealFreshWaterFlux ) sIceLoadFacLoc = recip_rhoConst
Expand Down Expand Up @@ -345,6 +360,28 @@ SUBROUTINE ECCO_PHYS( myTime, myIter, myThid )
ENDDO
areavolGlob=0. _d 0

#ifdef ALLOW_SEAICE
IF (gencost_barfile(kgen)(1:11).EQ.'m_freeboard') THEN
DO bj=myByLo(myThid),myByHi(myThid)
DO bi=myBxLo(myThid),myBxHi(myThid)
DO j = 1,sNy
DO i = 1,sNx
tmprecip_area = 1. _d 0/SQRT(area(i,j,bi,bj)*area(i,j,bi,bj)
& + area_reg_sq)
Copy link
Collaborator

Choose a reason for hiding this comment

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

because SEAICE_area_reg can be a read-in param from data.seaice, it's possible some user might set it to zero, thus resulting in division by zero here?

tmphsnowact = hsnow(I,J,bi,bj) * tmprecip_area
tmpheffact = heff(I,J,bi,bj) * tmprecip_area
gencost_storefld(i,j,bi,bj,kgen) =
& ( ( tmphsnowact + tmpheffact ) -
Copy link
Collaborator

Choose a reason for hiding this comment

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

could you add a comment somewhere to say what this equation says? it 's the full actual thickness of (ice+snow) minus the equivalent when converted to water?

& ( tmphsnowact * SEAICE_rhoSnow +
& tmpheffact * SEAICE_rhoIce ) * recip_rhoConst
& ) * maskC(i,j,1,bi,bj)
ENDDO
ENDDO
ENDDO
ENDDO
ENDIF
#endif

DO bj=myByLo(myThid),myByHi(myThid)
DO bi=myBxLo(myThid),myBxHi(myThid)
DO j = 1,sNy
Expand Down