C $Header: /u/gcmpack/MITgcm/pkg/seaice/seaice_init_fixed.F,v 1.30 2014/08/23 20:22:16 torge Exp $
C $Name: $
#include "SEAICE_OPTIONS.h"
CStartOfInterface
SUBROUTINE SEAICE_INIT_FIXED( myThid )
C *==========================================================*
C | SUBROUTINE SEAICE_INIT_FIXED
C | o Initialization of sea ice model.
C *==========================================================*
C *==========================================================*
IMPLICIT NONE
C === Global variables ===
#include "SIZE.h"
#include "EEPARAMS.h"
#include "PARAMS.h"
#include "GRID.h"
#include "FFIELDS.h"
#include "SEAICE_SIZE.h"
#include "SEAICE_PARAMS.h"
#include "SEAICE.h"
#include "SEAICE_TRACER.h"
C === Routine arguments ===
C myThid - Thread no. that called this routine.
INTEGER myThid
CEndOfInterface
C === Local variables ===
#ifdef SEAICE_ITD
C msgBuf :: Informational/error message buffer
CHARACTER*(MAX_LEN_MBUF) msgBuf
CHARACTER*15 HlimitMsgFormat
C k - loop counter for ITD categories
INTEGER k
_RL tmpVar
#endif
C i,j, bi,bj :: Loop counters
c INTEGER i, j, bi, bj
INTEGER kSurface
#ifdef ALLOW_SITRACER
INTEGER iTracer
#endif
#ifdef SHORTWAVE_HEATING
cif Helper variable for determining the fraction of sw radiation
cif penetrating the model shallowest layer
_RL dummyTime
_RL swfracba(2)
_RL tmpFac
#endif /* SHORTWAVE_HEATING */
IF ( buoyancyRelation .EQ. 'OCEANICP' ) THEN
kSurface = Nr
ELSE
kSurface = 1
ENDIF
C Initialize MNC variable information for SEAICE
IF ( useMNC .AND.
& (seaice_tave_mnc.OR.seaice_dump_mnc.OR.SEAICE_mon_mnc)
& ) THEN
CALL SEAICE_MNC_INIT( myThid )
ENDIF
C Only Master Thread updates parameter in commom block:
_BEGIN_MASTER(myThid)
C restart parameter
SEAICEmomStartBDF = 0
IF ( SEAICEuseBDF2 ) SEAICEmomStartBDF = nIter0
#ifdef SHORTWAVE_HEATING
tmpFac = -1.0
dummyTime = 1.0
swfracba(1) = ABS(rF(1))
swfracba(2) = ABS(rF(2))
CALL SWFRAC(
I 2, tmpFac,
U swfracba,
I dummyTime, 0, myThid )
SWFracB = swfracba(2)
#else /* SHORTWAVE_HEATING */
SWFracB = 0. _d 0
#endif /* SHORTWAVE_HEATING */
C-- Set mcPheePiston coeff (if still unset)
IF ( SEAICE_mcPheePiston.EQ.UNSET_RL ) THEN
IF ( SEAICE_availHeatFrac.NE.UNSET_RL ) THEN
SEAICE_mcPheePiston = SEAICE_availHeatFrac
& * drF(kSurface)/SEAICE_deltaTtherm
ELSE
SEAICE_mcPheePiston = MCPHEE_TAPER_FAC
& * STANTON_NUMBER * USTAR_BASE
SEAICE_mcPheePiston = MIN( SEAICE_mcPheePiston,
& drF(kSurface)/SEAICE_deltaTtherm )
ENDIF
ENDIF
C-- SItracer specifications for basic tracers
#ifdef ALLOW_SITRACER
DO iTracer = 1, SItrNumInUse
C "ice concentration" tracer that should remain .EQ.1.
IF (SItrName(iTracer).EQ.'one') THEN
SItrFromOcean0(iTracer) =ONE
SItrFromFlood0(iTracer) =ONE
SItrExpand0(iTracer) =ONE
SItrFromOceanFrac(iTracer) =ZERO
SItrFromFloodFrac(iTracer) =ZERO
ENDIF
C age tracer: no age in ocean, or effect from ice cover changes
IF (SItrName(iTracer).EQ.'age') THEN
SItrFromOcean0(iTracer) =ZERO
SItrFromFlood0(iTracer) =ZERO
SItrExpand0(iTracer) =ZERO
SItrFromOceanFrac(iTracer) =ZERO
SItrFromFloodFrac(iTracer) =ZERO
ENDIf
C salinity tracer:
IF (SItrName(iTracer).EQ.'salinity') THEN
SItrMate(iTracer) ='HEFF'
SItrExpand0(iTracer) =ZERO
IF ( SEAICE_salinityTracer ) THEN
SEAICE_salt0 = ZERO
SEAICE_saltFrac = ZERO
ENDIF
ENDIF
C simple, made up, ice surface roughness index prototype
IF (SItrName(iTracer).EQ.'ridge') THEN
SItrMate(iTracer) ='AREA'
SItrFromOcean0(iTracer) =ZERO
SItrFromFlood0(iTracer) =ZERO
SItrExpand0(iTracer) =ZERO
SItrFromOceanFrac(iTracer) =ZERO
SItrFromFloodFrac(iTracer) =ZERO
ENDIF
#ifdef SEAICE_GREASE
C grease ice tracer:
c (Smedrud and Martin, 2014, Ann. Glac.)
IF (SItrName(iTracer).EQ.'grease') THEN
SItrMate(iTracer) ='HEFF'
SItrFromOcean0(iTracer) =ZERO
SItrFromFlood0(iTracer) =ZERO
SItrExpand0(iTracer) =ZERO
SItrFromOceanFrac(iTracer) =ZERO
SItrFromFloodFrac(iTracer) =ZERO
ENDIF
#endif /* SEAICE_GREASE */
ENDDO
#endif /* ALLOW_SITRACER */
#ifdef SEAICE_ITD
C use Equ. 22 of Lipscomb et al. (2001, JGR) to generate ice
C thickness category limits:
C - dependends on given number of categories nITD
C - choose between original parameters of Lipscomb et al. (2001):
C c1=3.0/N, c2=15*c1, c3=3.0
C or emphasize thin end of ITD (in order to enhance ice growth):
C c1=1.5/N, c2=42*c1, c3=3.3
C -> HINT: set parameters c1, c2 and c3 in seaice_readparms.F
Hlimit(0) = 0. _d 0
IF ( nITD.GT.1 ) THEN
tmpVar = nITD
tmpVar = oneRL / tmpVar
Hlimit_c1 = Hlimit_c1*tmpVar
Hlimit_c2 = Hlimit_c2*Hlimit_c1
DO k=1,nITD-1
Hlimit(k) = Hlimit(k-1)
& + Hlimit_c1
& + Hlimit_c2
& *( oneRL + TANH( Hlimit_c3 *( FLOAT(k-1)*tmpVar - oneRL ) ) )
ENDDO
ENDIF
C thickest category is unlimited
Hlimit(nITD) = 999.9 _d 0
WRITE(msgBuf,'(A,I2,A)')
& ' SEAICE_INIT_FIXED: ', nITD,
& ' sea ice thickness categories'
CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
& SQUEEZE_RIGHT, myThid )
WRITE(HlimitMsgFormat,'(A,I2,A)') '(A,',nITD,'F6.2,F6.1)'
WRITE(msgBuf,HlimitMsgFormat)
& ' SEAICE_INIT_FIXED: Hlimit = ', (Hlimit(k),k=0,nITD)
CALL PRINT_MESSAGE( msgBuf, standardMessageUnit,
& SQUEEZE_RIGHT, myThid )
#endif /* SEAICE_ITD */
C Only Master Thread updates parameter in common block:
_END_MASTER(myThid)
#ifdef SEAICE_ALLOW_JFNK
C Only Master Thread updates parameter in common block:
_BEGIN_MASTER(myThid)
C initialise some diagnostic counters for the JFNK solver
totalNewtonIters = 0
totalNewtonFails = 0
totalKrylovIters = 0
totalKrylovFails = 0
totalJFNKtimeSteps = 0
_END_MASTER(myThid)
C this cannot be done here, because globalArea is only defined
C after S/R PACKAGES_INIT_FIXED, so we move it to S/R SEAICE_INIT_VARIA
CML CALL SEAICE_MAP2VEC( nVec, rAw, rAs,
CML & scalarProductMetric, .TRUE., myThid )
CML DO bj=myByLo(myThid),myByHi(myThid)
CML DO bi=myBxLo(myThid),myBxHi(myThid)
CML DO i=1,nVec
CML scalarProductMetric(i,1,bi,bj) =
CML & scalarProductMetric(i,1,bi,bj)/globalArea
CML ENDDO
CML ENDDO
CML ENDDO
#endif /* SEAICE_ALLOW_JFNK */
C-- all threads wait for master to finish initialisation of shared params
_BARRIER
#ifdef ALLOW_DIAGNOSTICS
IF ( useDiagnostics ) THEN
CALL SEAICE_DIAGNOSTICS_INIT( myThid )
ENDIF
#endif
C-- Summarise pkg/seaice configuration
CALL SEAICE_SUMMARY( myThid )
RETURN
END