Update of IS2 default parameters and associated change to regtest.

ukmo-waves · Jan 16, 2020 · f2a653b · f2a653b
1 parent deabe20
commit f2a653b
Show file tree

Hide file tree

Showing 8 changed files with 44 additions and 43 deletions.
diff --git a/model/ftn/w3sis2md.ftn b/model/ftn/w3sis2md.ftn
@@ -575,7 +575,7 @@ ELSE
       END SUBROUTINE INSIS2
 !/ ------------------------------------------------------------------- /
       SUBROUTINE W3SIS2 (A, DEPTH, CICE, ICEH, ICEF, ICEDMAX, IX, IY, &
-                         S, D, DCREEP, WN, CG, WN_R, CG_ICE, R)
+                         S, D, DISSIP, WN, CG, WN_R, CG_ICE, R)
 !/
 !/                  +-----------------------------------+
 !/                  | WAVEWATCH III           NOAA/NCEP |
@@ -619,7 +619,7 @@ ELSE
 !       IX,IY   Int   I     Not used
 !       S       R.A.  O     Source term (1-D version)
 !       D       R.A.  O     Diagonal part of scattering (1-D version)
-!       DCREEP  R.A.  O     Diagonal creep term (1-D version)
+!       DISSIP  R.A.  O     Diagonal dissipation term (1-D version)
 !       WN      R.A.  I     Wave number
 !       CG      R.A.  I     Group speed
 !       WN_R    R.A.  I     Wave number in ice
@@ -680,7 +680,7 @@ ELSE
 !/ Parameter list
       REAL, INTENT(IN)        :: A(NSPEC), DEPTH, CICE, ICEH, ICEDMAX
       INTEGER, INTENT(IN)     :: IX, IY
-      REAL, INTENT(OUT)       :: S(NSPEC), D(NSPEC), DCREEP(NSPEC), R(NK)
+      REAL, INTENT(OUT)       :: S(NSPEC), D(NSPEC), DISSIP(NSPEC), R(NK)
       REAL, INTENT(INOUT)     :: ICEF
       REAL, INTENT(IN)        :: WN(NK), CG(NK), WN_R(NK), CG_ICE(NK)
 
@@ -701,8 +701,6 @@ ELSE
       REAL, PARAMETER         :: YOUNG = 5.49E+9   ! Young modulus 
       REAL, PARAMETER         :: POISSON = 0.3     ! Poisson Ratio 
       REAL                    :: SIGMA_C
-!      REAL, PARAMETER         :: SIGMA_C =  2.7414E+05  ! Flexural strenght in Pa
-                                                         ! this is now a parametrer IS2PARS(19) 
       REAL, PARAMETER         :: DENS = 1025.0     ! ice density 
       REAL                    :: GAMMA_TOY 
       REAL, DIMENSION(NK)        :: WN_I, WN_RP, WSQ, WLG, WLG_I, CG_I,   &
@@ -717,7 +715,7 @@ ELSE
 !
       S = 0.
       D = 0.
-      DCREEP = 0.
+      DISSIP = 0.
       DUMMY = WN
       WN_I = 0.
       WN_RP = 0.
@@ -731,11 +729,9 @@ ELSE
       !B_COLE=1.205E-9 * EXP(0.55*1.60218E-19/(1.38064852e-23*263.15))
       B_COLE=1.205E-9 * EXP(IS2PARS(24)*1.60218E-19/(1.38064852e-23*268.15))
       TAU_D=B_COLE/0.07
-      S_D=LOG(SIG(:)*TAU_D)
-      !IF ((IX.EQ.50).AND.(IY.EQ.20)) WRITE(*,*) 'B', B_COLE, 'SD1' ,S_D(1), 'SD5', S_D(5), 'TAU_D', TAU_D
+      S_D=LOG(SIG(1:NK)*TAU_D)
       DELTA_D=IS2PARS(21) 
-      !SIGMA_RATIO =1.
-      !
+! 
       IF (IS2PARS(9).GT.0) ICEDMIN = IS2PARS(9)
       IF (IS2PARS(12).GT.0) THEN 
         B=IS2PARS(12)
@@ -951,29 +947,29 @@ ELSE
             END IF
           END IF  !end of test  (IS2PARS(3).GT.0.5)
 !
-!  5. creep dissipation 
+!  5. inelastic or anelastic dissipation 
 !      
         IF (IS2PARS(12).GT.0) THEN
           DO IK = 1, NK 
 !        
 ! The TANH((DMAX-D*WLG_I(IK))/DMAX*C) 
 ! is an ad hoc factor that goes to zero for WLG << DMAX and 1 for WLG >> DMAX 
-! this should probably be adjusted... 
+! this should probably be adjusted. 
 !
             IF (IS2PARS(14)*WLG_I(IK).LT.DMAX) THEN 
               SMOOTHD=TANH((DMAX-IS2PARS(14)*WLG_I(IK))/(DMAX*IS2PARS(13)))
-              IF (IS2PARS(23).LE.0.5) THEN
-                DCREEP(1+(IK-1)*NTH:IK*NTH)=CREEPFAC*4*CURV(IK)                & 
+              IF (IS2PARS(23).LE.0.5) THEN  ! this is the inelastic option 
+                DISSIP(1+(IK-1)*NTH:IK*NTH)=CREEPFAC*4*CURV(IK)                & 
                                             *((2*PI)/WLG_I(IK))**(IS2PARS(15)+1) &
                                             /(CGRATIO(IK)**1*R(IK)**2)           &
                                             *SMOOTHD
-              ELSE
-                DCREEP(1+(IK-1)*NTH:IK*NTH) =-4*4/3*SIG(IK)* DELTA_D*ALPHA_D *WN_I(IK)**4 * (YOUNG/MP2)**2         &
+              ELSE ! this is the inelastic option
+                DISSIP(1+(IK-1)*NTH:IK*NTH) =-4*4/3*SIG(IK)* DELTA_D*ALPHA_D *WN_I(IK)**4 * (YOUNG/MP2)**2         &
                                              * (ICEH/2)**3/3 * 1/( EXP(ALPHA_D*S_D(IK)) + EXP(-ALPHA_D*S_D(IK))) &
                                              * SMOOTHD /(R(IK)**2*CGRATIO(IK)) / (DWAT*GRAV) *TPIINV 
               END IF
             END IF
-            S=S+DCREEP*CICE*A
+            S=S+DISSIP*CICE*A
             END DO ! end of loop on IK
           ENDIF ! end of test (IS2PARS(12).GT.0)
 !

diff --git a/model/ftn/w3triamd.ftn b/model/ftn/w3triamd.ftn
@@ -221,6 +221,7 @@ CONTAINS
       ALLOCATE(XYBTMP1(NODES,3))
       DO I= 1, NODES
         READ(NDS,*) j, XYBTMP1(I,1), XYBTMP1(I,2), XYBTMP1(I,3)
+WRITE(6,*) NODES, I, j, XYBTMP1(I,1), XYBTMP1(I,2), XYBTMP1(I,3)
         END DO
 !
 ! read number of elements and elements from Gmsh files
@@ -730,6 +731,8 @@ CONTAINS
 !/ ------------------------------------------------------------------- /   
          USE W3GDATMD  
 !/S      USE W3SERVMD, ONLY: STRACE
+         USE W3ODATMD, ONLY: NDSE
+
          IMPLICIT NONE
 !	 
 !local parameters	 
@@ -764,7 +767,9 @@ CONTAINS
          I2 = TRIGP(K,2)
          I3 = TRIGP(K,3)
          TRIA(K) = -1.d0*TRIA(K)
-         STOP 'WRONG TRIANGLE' 
+         WRITE(NDSE,*) 'WRONG TRIANGLE',TRIA(K),K,I1,I2,I3, XYB(I2,2)-XYB(I1,2), &
+                          XYB(I1,1)-XYB(I3,1),XYB(I3,2)-XYB(I1,2), XYB(I2,1)-XYB(I1,1)
+         STOP 
          END IF    
        END DO
      END SUBROUTINE

diff --git a/model/ftn/ww3_grid.ftn b/model/ftn/ww3_grid.ftn
@@ -2210,11 +2210,11 @@
       I1E3D=1
       I2E3D=NK
       P2SF   = 0
-      I1P2SF = 1
+      I1P2SF = 1	
       I2P2SF = 15
-      US3D=0 
-      I1US3D=1
-      I2US3D=NK
+      US3D   = 0 
+      I1US3D = 1
+      I2US3D = NK
       USSP=0
       IUSSP=1
       STK_WN(:)=0.0
@@ -2255,15 +2255,15 @@
 !/IS2      IS2FRAGILITY=0.9
 !/IS2      IS2DMIN=20
 !/IS2      IS2DAMP=0.
-!/IS2      IS2CREEPB=0.
-!/IS2      IS2CREEPC=0.2
-!/IS2      IS2CREEPD=0.4
+!/IS2      IS2CREEPB=2.0E+07
+!/IS2      IS2CREEPC=0.4
+!/IS2      IS2CREEPD=0.5
 !/IS2      IS2CREEPN=3.0
 !/IS2      IS2BREAKE=1.
 !/IS2      IS2WIM1=1.
-!/IS2      IS2ANDISB=.FALSE. !anelastic instead of inelastic dissipation if IS2CREEPB>0
+!/IS2      IS2ANDISB=.TRUE.  !anelastic instead of inelastic dissipation if IS2CREEPB>0
 !/IS2      IS2ANDISE=0.55    !energy of activation
-!/IS2      IS2ANDISD=1.672E-9!
+!/IS2      IS2ANDISD=2.0E-9  !relaxation of dislocations in 1/Pa
 !/IS2      IS2ANDISN=1.      !dependency on stress. Equal to 1 normally?
 !/IS2      CALL READNL ( NDSS, 'SIS2', STATUS )
 !/IS2      WRITE (NDSO,947) STATUS

diff --git a/regtests/ww3_tic1.4/input/ww3_ounf.inp b/regtests/ww3_tic1.4/input/ww3_ounf.inp
@@ -4,7 +4,7 @@ $--------------------------------------------------------------------- $
 $ First output time (yyyymmdd hhmmss), increment of output (s), 
 $ and number of output times.
 $
-  19680606 000000  3600.  1000
+  20120101 000000  3600.  1000
 $
 $ Fields requested --------------------------------------------------- $
 $

diff --git a/regtests/ww3_tic1.4/input/ww3_shel.inp b/regtests/ww3_tic1.4/input/ww3_shel.inp
@@ -1,10 +1,10 @@
 $ WAVEWATCH III shell input file
 $ ------------------------------
-   T F     Ice parameter 1
+   T T     Ice parameter 1
    T T     Ice parameter 2
    F F     Ice parameter 3
    F F     Ice parameter 4
-   T F     Ice parameter 5
+   F F     Ice parameter 5
    F F     Mud parameter 1
    F F     Mud parameter 2
    F F     Mud parameter 3
@@ -16,12 +16,12 @@ $ ------------------------------
    F       Assimilation data : 1-D spectra
    F       Assimilation data : 2-D spectra.
 $
-   20120101 000000
-   20120102 120000
+   20120103 180000
+   20120105 120000
 $
    1
 $
-   20120101 000000    900  20120102 120000
+   20120101 000000    900  20120108 120000
 N
 $ Options : DPT CUR WND AST WLV ICE IBG D50 IC1 IC5 HS LM T02 T0M1 T01 FP DIR SPR
 $           DP HIG EF TH1M STH1M TH2M STH2M WN PHS PTP PLP PDIR PSPR PWS PDP
@@ -64,7 +64,7 @@ $ Homogeneous field data --------------------------------------------- $
 $      ICE1  ice thickness, e.g. 1 m
 $      ICE2  eddy viscosity parameter of Liu et al.
 $ constant case:
-$   'IC1' 20120101 000000   0.1
+   'IC1' 20120101 000000   0.8
    'IC2' 20120101 000000   1536.0e-4
 $   'IC5' 20120101 000000   30.
    'LEV' 20120101 000000   1.1

diff --git a/regtests/ww3_tic1.4/input/ww3_shel.nml b/regtests/ww3_tic1.4/input/ww3_shel.nml
@@ -58,9 +58,9 @@
 !     INPUT%ASSIM%SPEC2D          = 'F'
 ! -------------------------------------------------------------------- !
 &INPUT_NML
-  INPUT%FORCING%ICE_PARAM1    = 'T'
-  INPUT%FORCING%ICE_PARAM2    = 'H'
-  INPUT%FORCING%ICE_PARAM5    = 'T'
+  INPUT%FORCING%ICE_PARAM1    = 'H'
+  INPUT%FORCING%ICE_PARAM2    = 'F'
+  INPUT%FORCING%ICE_PARAM5    = 'F'
   INPUT%FORCING%WATER_LEVELS  = 'H'
   INPUT%FORCING%ICE_CONC      = 'T'
 /
@@ -330,14 +330,14 @@
 !     HOMOG_INPUT(I)%VALUE3         =  0
 ! -------------------------------------------------------------------- !
 &HOMOG_COUNT_NML
-  HOMOG_COUNT%N_IC2                =  1
+  HOMOG_COUNT%N_IC1                =  1
   HOMOG_COUNT%N_LEV                =  1
 /
 
 &HOMOG_INPUT_NML
-  HOMOG_INPUT(1)%NAME        = 'IC2'
+  HOMOG_INPUT(1)%NAME        = 'IC1'
   HOMOG_INPUT(1)%DATE        = '20120101 000000'
-  HOMOG_INPUT(1)%VALUE1      = 1536.0e-4
+  HOMOG_INPUT(1)%VALUE1      = 1.
 
   HOMOG_INPUT(2)%NAME        = 'LEV'
   HOMOG_INPUT(2)%DATE        = '20120101 000000'

diff --git a/regtests/ww3_tic1.4/input/ww3_strt.inp b/regtests/ww3_tic1.4/input/ww3_strt.inp
@@ -5,8 +5,8 @@ $  0.1  0.0001 225. 12   0.    -5.E3  0.  5.E3  1.0
 $  0.1  0.0001 315. 12   0.    -5.E3  0.  5.E3  1.0
 $  0.1  0.0001 240.  2   0.    -5.E3  0.  5.E3  1.0
 $  fp   sip    thm  ncos xm     six   ym  siy  hmax
-$   0.07  0.05   270. 100   0.    -5.E3  0.  5.E3  3.0
-   0.1  0.03   270. 100   0.    -5.E3  0.  5.E3  3.0
+$   0.07  0.05   270. 100   0.    -5.E3  0.  5.E3  2.0
+   0.1  0.02   270. 100   0.    -5.E3  0.  5.E3  3.0
 $  0.1  0.0001 270. 100  0.    -5.E3  0.  5.E3  3.0
 $  0.1  0.0001 270.  2   0.    -5.E3  0.  5.E3  1.0
 $  0.1  0.0001 300.  2   0.    -5.E3  0.  5.E3  1.0

diff --git a/regtests/ww3_tic2.1/info b/regtests/ww3_tic2.1/info
@@ -46,7 +46,7 @@
 #  Example run_test command(s):                                             #
 #  (some details will vary by local system and configuration)               #
 #  ./bin/run_test -w work_IC1 -i input_IC1 ../model ww3_tic2.1              #
-#  ./bin/run_test -w work_IC2IS2 -c gnu_debug -i input_IC2IS2               #
+#  ./bin/run_test -w work_IC2IS2 -c gnu_debug -i input_IC2IS2  -o netcdf    #
 #                                     ../model ww3_tic2.1                   #
 #  ./bin/run_test -p mpirun -n 3 -w work_IC4 -i input_IC4                   #
 #                                     ../model ww3_tic2.1                   #