!> @file !> @brief Floe-size dependant scattering of waves in the marginal ice zone. !> @details Based on tabulated scattering coefficients for a semi-infinite !> ice sheet. !> !> @author F. Ardhuin !> @author P. Nicot !> @author C. Sevigny !> @author G. Boutin !> @date 21-Jan-2018 !> #include "w3macros.h" !/ ------------------------------------------------------------------- / !> !> @brief Floe-size dependant scattering of waves in the marginal ice zone. !> !> @details based on tabulated scattering coefficients for a semi-infinite !> ice sheet. See papers by Dumont et al. (JGR 2011) and Williams et al. !> (OM 2013) combined with flexural dissipation and ice break-up. !> !> @author F. Ardhuin !> @author P. Nicot !> @author C. Sevigny !> @author G. Boutin !> @date 21-Jan-2018 !> !> @copyright Copyright 2009-2022 National Weather Service (NWS), !> National Oceanic and Atmospheric Administration. All rights !> reserved. WAVEWATCH III is a trademark of the NWS. !> No unauthorized use without permission. !> MODULE W3SIS2MD !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | F. Ardhuin | !/ | P. Nicot | !/ | C. Sevigny | !/ | G. Boutin | !/ | FORTRAN 90 | !/ | Last update : 21-Jan-2018 | !/ +-----------------------------------+ !/ !/ For updates see W3SID1 documentation. !/ ! 1. Purpose : ! ! Floe-size dependant scattering of waves in the marginal ice zone based on tabulated ! scattering coefficients for a semi-infinite ice sheet. See papers ! by Dumont et al. (JGR 2011) and Williams et al. (OM 2013) ! combined with flexural dissipation and ice break-up. ! ! 2. Variables and types : ! ! 3. Subroutines and functions : ! ! Name Type Scope Description ! ---------------------------------------------------------------- ! W3SIS2 Subr. Public Ice scattering term. ! ---------------------------------------------------------------- ! ! 4. Subroutines and functions used : ! ! See subroutine documentation. ! ! 5. Remarks : ! ! 6. Switches : ! ! See subroutine documentation. ! ! 7. Source code : !/ !/ ------------------------------------------------------------------- / !/ INTEGER , PARAMETER :: NTHICK = 20, NICED = 500 REAL, PARAMETER :: FRAGILITY = 0.9 REAL :: THICK1 = 0.1, DTHICK = 0.25 REAL :: ICEDMIN ! minimum floe diameter REAL :: ICEDAVETAB(NICED) REAL, DIMENSION(:) , ALLOCATABLE :: SIS2ALPHAS(:,:),SIS2ALPHA2(:,:) DOUBLE PRECISION, ALLOCATABLE,DIMENSION(:,:) :: IS2EIGVEC, IS2SCATMAT DOUBLE PRECISION , ALLOCATABLE,DIMENSION(:) :: IS2EIGVAL PRIVATE :: SIS2ALPHAS, SIS2ALPHA2 PUBLIC :: IS2EIGVEC, IS2EIGVAL ! PUBLIC :: W3SIS2, INSIS2, W3RPWNICE PRIVATE :: FINDROOTS_NR, W3FSD_DAVE, FUNCD_FVAL, FUNCD_FDERIV !/ CONTAINS !/ ------------------------------------------------------------------- / !> !> @brief Fill tables used for scattering. !> !> @details Linear interpolation. !> !> @author P. Nicot !> @author F. Ardhuin !> @date 21-Jan-2018 !> SUBROUTINE INSIS2 !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | P. Nicot & F. Ardhuin | !/ | FORTRAN 90 | !/ | Last update : 21-Jan-2018 | !/ +-----------------------------------+ !/ !/ 01-Apr-2014 : Creation ( version 4.18 ) !/ 13-Dec-2015 : Adds diagonalization of scat. matrix( version 5.10 ) !/ 21-Jan-2018 : Implements non-isotropic example ( version 6.04 ) !/ ! 1. Purpose : ! ! Fill tables used for scattering ! ! 2. Method : ! ! Linear interpolation ! ! 3. Parameters : ! ! See module documentation. ! ! 4. Error messages : ! ! - None. ! ! 5. Called by : ! ! - W3IOGR (initialization after reading mod_def.ww3) ! ! 6. Subroutines used : ! ! - None ! ! 7. Remarks : ! ! ! 8. Structure ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / !/ USE W3GDATMD, ONLY: SIG, DSIP, NK, NTH, IS2PARS, & EC2, ES2, ESC, ESIN, ECOS USE CONSTANTS, ONLY: TPI, TPIINV USE W3SERVMD, ONLY: DIAGONALIZE #ifdef W3_S USE W3SERVMD, ONLY: STRACE #endif ! IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ INTEGER :: I, J, K, IND, NFTAB, NROT REAL :: SIS1HTABLE(20), SIS1FTABLE(25) REAL :: SIS1ALPHATABLE(NTHICK,25), X REAL :: SIS1ALPHATABLE2(NTHICK,25) #ifdef W3_S INTEGER, SAVE :: IENT = 0 #endif !/ !/ ------------------------------------------------------------------- / !/ #ifdef W3_S CALL STRACE (IENT, 'SIS2ALPHATAB') #endif ! ! -------------------------------------------------------------------- / ! 1. Fills array of reflection as a function of frequency and ice thickness ! ALLOCATE(SIS2ALPHAS(NTHICK,NK)) ALLOCATE(SIS2ALPHA2(NTHICK,NK)) ! ! Table of ice thickness for which the reflection was computed ! SIS1HTABLE = (/ 0.1, 0.35, 0.6, 0.85, 1.1, 1.35, 1.6, 1.85, 2.1, 2.35, & 2.6, 2.85, 3.1, 3.35, 3.6, 3.85, 4.1, 4.35, 4.6, 4.85 /) NFTAB = 25 ! ! Table of frequencies for which the reflection was computed ! SIS1FTABLE = (/ 0.0420, 0.04620, 0.050820, 0.0559020, 0.06149220, 0.06764142,0.0744055620000000, & 0.0818461182, 0.09003073002, 0.099033803022, 0.1089371833242, 0.11983090, 0.13181399, & 0.144995391, 0.15949493 , 0.175444423115458, 0.192988865427003, 0.212287751969704, & 0.233516527, 0.256868, 0.28255499787, 0.310810497658843, 0.341891547424728, 0.376080702167200, 0.413688772383920 /) IF (IS2PARS(18).LT.0.5) THEN SIS1ALPHATABLE = reshape((/ & 1.78E-007, 2.21E-006, 6.57E-006, 1.31E-005, 2.28E-005, 3.60E-005, 5.23E-005, 7.19E-005, 9.60E-005, 0.0001260665, & 0.0001621645, 0.0002032202, 0.0002483457, 0.0002987263, 0.0003571903, & 0.0004247358, 0.00049714, 0.0005689017, 0.0006469729, 0.0007470985, & 2.74E-007, 3.40E-006, 1.04E-005, 2.10E-005, 3.67E-005, 5.82E-005, & 8.52E-005, 0.0001179843, 0.0001581719, 0.0002075378, & 0.0002663323, 0.0003330095, 0.0004063158, 0.000487845, 0.0005815828, & 0.0006893214, 0.0008061048, 0.0009252626, 0.0010565619, 0.0012218782, & 4.16E-007, 5.21E-006, 1.62E-005, 3.32E-005, 5.91E-005, 9.56E-005, 0.0001423817, & 0.0002001783, 0.0002716622, 0.0003599526, 0.0004656353, & 0.0005861828, 0.0007194785, 0.0008683691, 0.0010401368, 0.0012386623, 0.0014562639, & 0.0016819935, 0.0019335276, 0.0022507523, & 6.25E-007, 7.95E-006, 2.50E-005, 5.28E-005, 9.71E-005, 0.0001615969, 0.0002467311, & 0.0003540015, 0.0004891301, 0.0006588358, 0.0008648602, & 0.0011024698, 0.0013675698, 0.0016664349, 0.0020150216, 0.0024225002, 0.0028732348, & 0.0033437743, 0.0038725363, 0.0045480433, & 9.32E-007, 1.21E-005, 3.91E-005, 8.60E-005, 0.000165062, 0.0002848805, 0.0004473356, & 0.0006565411, 0.0009257279, 0.0012707025, 0.0016966237, & 0.0021939892, 0.002754291, 0.0033926056, 0.0041466541, 0.0050385513, 0.0060311467, & 0.0070669538, 0.0082363302, 0.0097527074, & 1.38E-006, 1.86E-005, 6.30E-005, 0.0001464051, 0.0002950074, 0.0005287398, 0.0008535728, & 0.0012801516, 0.0018391215, 0.0025674934, 0.0034790152,& 0.0045539812, 0.0057742581, 0.0071757314, 0.0088467499, 0.0108382999, & 0.0130588683, 0.0153657347, 0.0179697664, 0.0213785773, & 2.05E-006, 2.91E-005, 0.0001064725, 0.0002646727, 0.0005601397, 0.0010380532, & 0.0017148531, 0.0026159282, 0.0038087163, 0.0053744112, & 0.0073439386, 0.0096752268, 0.0123295445, 0.0153861547, 0.0190373907, 0.0233881203, & 0.0282220195, 0.0332084504, 0.0388053869, 0.0461221517, & 3.03E-006, 4.71E-005, 0.0001918323, 0.000513363, 0.0011329164, & 0.00214921, 0.0036015483, 0.0055437802, 0.0081104965, 0.0114597841, 0.0156453149, & 0.0205751599, 0.0261681604, 0.0325784506, 0.040175062, 0.0491364185, & 0.0589976693, 0.0690831851, 0.0802738154, 0.0946821772, & 4.53E-006, 8.02E-005, 0.0003720589, 0.001067721, 0.0024215478, & 0.0046335423, 0.0077820302, 0.0119609444, 0.0173982353, 0.0243482166, & 0.0328666235, 0.0427569128, 0.0538584883, 0.0664227317, 0.0810452442, & 0.0979544205, 0.1162816559, 0.1348568677, 0.1551243874, 0.1804646861, & 6.84E-006, 0.0001465372, 0.0007768361, 0.0023492877, 0.0053624751, & 0.0101694228, 0.0168827981, 0.0256087065, 0.0366354159, 0.0502611381, & 0.0664627541, 0.0848577616, 0.1051596929, 0.1276991046, 0.1532607056, & 0.1820303259, 0.2126339282, 0.2433779188, 0.276250064, 0.3157667291, & 1.06E-005, 0.0002908265, 0.0017246725, 0.0053239059, 0.0119340853, & 0.0220417832, 0.035687347, 0.0528302869, 0.073649217, 0.098327901, & 0.1266340957, 0.1579155214, 0.1917144685, 0.2283923343, 0.2688390962, & 0.3131039168, 0.359282913, 0.4052034867, 0.4533028057, 0.5090275342, & 1.68E-005, 0.0006285069, 0.0039687874, 0.0119764661, 0.0257390989, & 0.0456577018, 0.0713257381, 0.1021551099, 0.1379813612, 0.1787934975, & 0.2240723105, 0.2727534775, 0.3241066558, 0.3786049119, 0.4374575502, & 0.5006384418, 0.5653760852, 0.6286093824, 0.6937583927, 0.7684424911, & 2.82E-005, 0.0014590465, 0.0091291531, 0.0257006983, 0.051998859, & 0.0879406127, 0.1317806874, 0.181756664, 0.2375759087, 0.2995488732, & 0.3668698493, 0.4374715634, 0.5099661328, 0.5855233305, 0.6668652222, & 0.754201342, 0.8422038741, 0.925009729, 1.0096608603, 1.1112763605, & 5.06E-005, 0.0035218608, 0.0201532593, 0.0509283911, 0.0962607638, & 0.1553638291, 0.2233707746, 0.2967358894, 0.3765436929, 0.4650074273, & 0.5607827923, 0.6590228472, 0.756685343, 0.8572097524, 0.9677298862,& 1.0892101997, 1.2094355521, 1.3151711317, 1.4236588573, 1.5698199322, & 9.99E-005, 0.0084186275, 0.0412091805, 0.0918067318, 0.1629311719, & 0.2530998556, 0.3511295005, 0.4513397894, 0.5593550205, 0.6817564883, & 0.815653431, 0.9499632903, 1.0782538468, 1.2090808019, 1.3588213156, & 1.5298071338, 1.6953922745, 1.8273584757, 1.9646959304, 2.1824935098, & 0.0002197702, 0.0188788997, 0.076522568, 0.1517352519, 0.2567566667, & 0.3885311336, 0.5247441812, 0.6573850155, 0.8007012438, 0.9681737964, & 1.1534672992, 1.3345908371, 1.4996817772, 1.666222509, 1.8655897771, & 2.1018122262, 2.3245234997, 2.4819251275, 2.6492217053, 2.9647668115, & 0.0005376597, 0.0381101752, 0.1293683129, 0.2359311434, 0.3874916537, & 0.576671054, 0.7645781239, 0.9404856491, 1.1305780655, 1.356723506, & 1.6068566622, 1.844140891, 2.0502598744, 2.2550088512, 2.5085040228, & 2.815733869, 3.096422752, 3.2715199709, 3.462031842, 3.8799520309, & 0.0014123358, 0.0685490628, 0.203577783, 0.354894301, 0.5714749709, & 0.8380065625, 1.0957796855, 1.3307089678, 1.5810360789, 1.8762692254,& 2.197628874, 2.4940515596, 2.7425260367, 2.9852393578, 3.2886471617, & 3.6568112069, 3.9833963075, 4.1693103781, 4.3750188868, 4.868451483, & 0.0037199175, 0.1126937266, 0.3066036571, 0.5225755275, 0.8226663225, & 1.1819668289, 1.5247740065, 1.8332140124, 2.152798111, 2.5177593097, & 2.9057350396, 3.2594397982, 3.5553142762, 3.8420910788, 4.1929464318,& 4.6097641365, 4.975238806, 5.1850270882, 5.4153321104, 5.9516719058, & 0.009037082, 0.1761561171, 0.4486599083, 0.747307698, 1.1406376531, & 1.5966776005, 2.0299434584, 2.4186982894, 2.8096987013, 3.2399866752, & 3.6899433755, 4.105747098, 4.4658248514, 4.8172532951, 5.2300901821, & 5.7046450426, 6.1270633187, 6.3975380116, 6.6851833028, 7.2647998139, & 0.0191118494, 0.2660897811, 0.6359174284, 1.0316197503, 1.5313672569, & 2.0961171527, 2.6294074335, 3.1054178673, 3.5736395046, 4.0762156975, & 4.5979660806, 5.0883063988, 5.5262800229, 5.9567260476, 6.4472029911, & 6.9971641707, 7.4947435133, 7.841785986, 8.2023388524, 8.8466621683, & 0.0352335589, 0.3840118688, 0.8727819679, 1.3979077636, 2.0496249898, & 2.7683726012, 3.4325758193, 4.0133546648, 4.572769913, 5.1624475081, & 5.7663362962, 6.3283500032, 6.8264493734, 7.3114532863, 7.8576036504, & 8.4626783134, 9.0041860776, 9.378581111, 9.7670030485, 10.46036883, & 0.0588123086, 0.5407217302, 1.1866725247, 1.8744101474, 2.7029832853, & 3.586993579, 4.3802119773, 5.0573372156, 5.6991092734, 6.3694458919, & 7.0521552001, 7.6857035561, 8.2469182713, 8.7926955021, 9.4042992778, & 10.0780089999, 10.6798375161, 11.0993158393, 11.5359791079, 12.3054058274, & 0.0912925253, 0.7587053794, 1.5933132286, 2.4467380656, 3.4467575836, & 4.4845019618, 5.3910117536, 6.1520041613, 6.8749329856, 7.6395860278, & 8.4275368684, 9.1660844666, 9.8268448897, 10.4746437466, 11.2031746957, & 12.0078539018, 12.7326579385, 13.2484048884, 13.7868789237, 14.7164075718, & 0.1320035456, 1.258347597, 2.2697363962, 2.8626691529, 3.4800285532, & 4.1591450976, 4.8036627318, 5.4044338329, 6.0517560386, 6.8415021436, & 7.8080659794, 8.9257124242, 10.1640694041, 11.5476843252, 13.1494718506, & 14.9965430946, 16.972215146, 18.903468102, 20.9125176371, 23.3776351255 & /) ,(/NTHICK,NFTAB/)) ELSE ! May be changed, but according to T. Williams, wim1 is okay from 0.25 only SIS1HTABLE = (/ 0.25, 0.35, 0.6, 0.85, 1.1, 1.35, 1.6, 1.85, 2.1, 2.35, & 2.6, 2.85, 3.1, 3.35, 3.6, 3.85, 4.1, 4.35, 4.6, 4.85 /) SIS1ALPHATABLE = reshape((/ & 3.80373e-06 , 6.02822e-06 , 1.12121e-05 , 2.24588e-05 , & 3.05165e-05 , 3.85142e-05 , 5.29712e-05 , 7.55453e-05 , & 1.05531e-04 , 1.38952e-04 , 1.73378e-04 , 2.09421e-04 , & 2.47299e-04 , 2.87669e-04 , 3.34304e-04 , 3.91094e-04 , & 4.57710e-04 , 5.31901e-04 , 6.14490e-04 , 7.08299e-04 , & 4.76980e-06 , 6.95433e-06 , 1.54179e-05 , 2.46138e-05 , & 3.79955e-05 , 5.96725e-05 , 8.81668e-05 , 1.20881e-04 , & 1.57364e-04 , 1.99666e-04 , 2.50723e-04 , 3.13346e-04 , & 3.89907e-04 , 4.81606e-04 , 5.88164e-04 , 7.08707e-04 , & 8.42536e-04 , 9.88622e-04 , 1.14498e-03 , 1.30931e-03 , & 5.84556e-06 , 8.58738e-06 , 2.00947e-05 , 3.18410e-05 , & 5.67385e-05 , 9.37456e-05 , 1.37304e-04 , 1.87832e-04 , & 2.49218e-04 , 3.28900e-04 , 4.31988e-04 , 5.59364e-04 , & 7.11302e-04 , 8.86736e-04 , 1.08074e-03 , 1.28808e-03 , & 1.50809e-03 , 1.74230e-03 , 1.98930e-03 , 2.24607e-03 , & 6.90305e-06 , 1.12155e-05 , 2.59514e-05 , 5.15063e-05 , & 9.18902e-05 , 1.45703e-04 , 2.18304e-04 , 3.18476e-04 , & 4.51992e-04 , 6.21037e-04 , 8.24428e-04 , 1.05881e-03 , & 1.32112e-03 , 1.60970e-03 , 1.92396e-03 , 2.26484e-03 , & 2.63577e-03 , 3.04238e-03 , 3.49116e-03 , 3.98911e-03 , & 8.46004e-06 , 1.46791e-05 , 3.69837e-05 , 8.32870e-05 , & 1.49151e-04 , 2.46175e-04 , 3.92165e-04 , 5.95871e-04 , & 8.57938e-04 , 1.17280e-03 , 1.53668e-03 , 1.95109e-03 , & 2.42048e-03 , 2.95341e-03 , 3.56468e-03 , 4.27150e-03 , & 5.08843e-03 , 6.02927e-03 , 7.11153e-03 , 8.35459e-03 , & 1.11513e-05 , 1.95090e-05 , 5.78215e-05 , 1.32267e-04 , & 2.57932e-04 , 4.61051e-04 , 7.55044e-04 , 1.14193e-03 , & 1.62478e-03 , 2.21208e-03 , 2.92101e-03 , 3.77684e-03 , & 4.80845e-03 , 6.04635e-03 , 7.52256e-03 , 9.26818e-03 , & 1.13101e-02 , 1.36712e-02 , 1.63732e-02 , 1.94360e-02 , & 1.52571e-05 , 2.78128e-05 , 9.43665e-05 , 2.30636e-04 , & 4.90057e-04 , 9.02105e-04 , 1.48012e-03 , 2.24822e-03 , & 3.24754e-03 , 4.53444e-03 , 6.16979e-03 , 8.21171e-03 , & 1.07145e-02 , 1.37253e-02 , 1.72805e-02 , 2.14091e-02 , & 2.61367e-02 , 3.14840e-02 , 3.74632e-02 , 4.40798e-02 , & 2.13731e-05 , 4.29771e-05 , 1.64050e-04 , 4.49510e-04 , & 9.84576e-04 , 1.82921e-03 , 3.06425e-03 , 4.79802e-03 , & 7.14520e-03 , 1.02139e-02 , 1.40940e-02 , 1.88529e-02 , & 2.45396e-02 , 3.11875e-02 , 3.88141e-02 , 4.74255e-02 , & 5.70233e-02 , 6.76049e-02 , 7.91611e-02 , 9.16776e-02 , & 3.17688e-05 , 7.05749e-05 , 3.13449e-04 , 9.26652e-04 , & 2.08059e-03 , 3.98338e-03 , 6.87070e-03 , 1.09528e-02 , & 1.63868e-02 , 2.32726e-02 , 3.16649e-02 , 4.15851e-02 , & 5.30293e-02 , 6.59789e-02 , 8.04085e-02 , 9.62863e-02 , & 1.13575e-01 , 1.32236e-01 , 1.52233e-01 , 1.73533e-01 , & 5.14019e-05 , 1.24975e-04 , 6.49063e-04 , 2.01553e-03 , & 4.70028e-03 , 9.17913e-03 , 1.58010e-02 , 2.47576e-02 , & 3.61186e-02 , 4.98708e-02 , 6.59568e-02 , 8.42994e-02 , & 1.04811e-01 , 1.27404e-01 , 1.51996e-01 , 1.78511e-01 , & 2.06873e-01 , 2.37011e-01 , 2.68862e-01 , 3.02365e-01 , & 9.05394e-05 , 2.44433e-04 , 1.43884e-03 , 4.66981e-03 , & 1.09558e-02 , 2.08947e-02 , 3.46681e-02 , 5.22109e-02 , & 7.33489e-02 , 9.78738e-02 , 1.25579e-01 , 1.56274e-01 , & 1.89790e-01 , 2.25973e-01 , 2.64687e-01 , 3.05808e-01 , & 3.49223e-01 , 3.94830e-01 , 4.42530e-01 , 4.92227e-01 , & 1.75074e-04 , 5.25080e-04 , 3.39573e-03 , 1.10036e-02 , & 2.46312e-02 , 4.44096e-02 , 6.99511e-02 , 1.00742e-01 , & 1.36302e-01 , 1.76219e-01 , 2.20147e-01 , 2.67802e-01 , & 3.18937e-01 , 3.73339e-01 , 4.30813e-01 , 4.91183e-01 , & 5.54282e-01 , 6.19954e-01 , 6.88048e-01 , 7.58418e-01 , & 3.74908e-04 , 1.22443e-03 , 8.17252e-03 , 2.46994e-02 , & 5.11207e-02 , 8.63422e-02 , 1.29173e-01 , 1.78641e-01 , & 2.34006e-01 , 2.94691e-01 , 3.60235e-01 , 4.30258e-01 , & 5.04426e-01 , 5.82438e-01 , 6.64019e-01 , 7.48914e-01 , & 8.36872e-01 , 9.27649e-01 , 1.02102e+00 , 1.11676e+00 , & 8.82800e-04 , 3.03043e-03 , 1.88248e-02 , 5.09125e-02 , & 9.68861e-02 , 1.54064e-01 , 2.20524e-01 , 2.94949e-01 , & 3.76394e-01 , 4.64138e-01 , 5.57588e-01 , 6.56231e-01 , & 7.59605e-01 , 8.67281e-01 , 9.78849e-01 , 1.09393e+00 , & 1.21214e+00 , 1.33315e+00 , 1.45662e+00 , 1.58224e+00 , & 2.22521e-03 , 7.51676e-03 , 3.98357e-02 , 9.57376e-02 , & 1.68983e-01 , 2.55527e-01 , 3.52939e-01 , 4.59644e-01 , & 5.74490e-01 , 6.96566e-01 , 8.25086e-01 , 9.59331e-01 , & 1.09864e+00 , 1.24238e+00 , 1.38997e+00 , 1.54085e+00 , & 1.69451e+00 , 1.85049e+00 , 2.00833e+00 , 2.16765e+00 , & 5.66833e-03 , 1.75464e-02 , 7.66453e-02 , 1.65727e-01 , & 2.74951e-01 , 3.99380e-01 , 5.36249e-01 , 6.83721e-01 , & 8.40370e-01 , 1.00497e+00 , 1.17642e+00 , 1.35368e+00 , & 1.53580e+00 , 1.72190e+00 , 1.91118e+00 , 2.10291e+00 , & 2.29645e+00 , 2.49123e+00 , 2.68675e+00 , 2.88259e+00 , & 1.36539e-02 , 3.72456e-02 , 1.35140e-01 , 2.67780e-01 , & 4.22901e-01 , 5.95111e-01 , 7.81333e-01 , 9.79307e-01 , & 1.18712e+00 , 1.40304e+00 , 1.62548e+00 , 1.85300e+00 , & 2.08433e+00 , 2.31832e+00 , 2.55400e+00 , 2.79055e+00 , & 3.02727e+00 , 3.26361e+00 , 3.49910e+00 , 3.73337e+00 , & 2.98722e-02 , 7.15314e-02 , 2.21423e-01 , 4.09271e-01 , & 6.21810e-01 , 8.53369e-01 , 1.10028e+00 , 1.35949e+00 , & 1.62829e+00 , 1.90425e+00 , 2.18522e+00 , 2.46937e+00 , & 2.75516e+00 , 3.04136e+00 , 3.32695e+00 , 3.61118e+00 , & 3.89345e+00 , 4.17331e+00 , 4.45047e+00 , 4.72471e+00 , & 5.88578e-02 , 1.25612e-01 , 3.41821e-01 , 5.98292e-01 , & 8.81793e-01 , 1.18598e+00 , 1.50596e+00 , 1.83748e+00 , & 2.17675e+00 , 2.52062e+00 , 2.86649e+00 , 3.21233e+00 , & 3.55662e+00 , 3.89823e+00 , 4.23636e+00 , 4.57049e+00 , & 4.90029e+00 , 5.22559e+00 , 5.54632e+00 , 5.86250e+00 , & 1.05468e-01 , 2.04724e-01 , 5.03158e-01 , 8.44026e-01 , & 1.21416e+00 , 1.60557e+00 , 2.01139e+00 , 2.42573e+00 , & 2.84383e+00 , 3.26205e+00 , 3.67776e+00 , 4.08914e+00 , & 4.49503e+00 , 4.89472e+00 , 5.28788e+00 , 5.67440e+00 , & 6.05432e+00 , 6.42782e+00 , 6.79514e+00 , 7.15656e+00 , & 1.74545e-01 , 3.14137e-01 , 7.13063e-01 , 1.15691e+00 , & 1.63126e+00 , 2.12499e+00 , 2.62866e+00 , 3.13493e+00 , & 3.63860e+00 , 4.13623e+00 , 4.62576e+00 , 5.10610e+00 , & 5.57680e+00 , 6.03787e+00 , 6.48955e+00 , 6.93227e+00 , & 7.36652e+00 , 7.79283e+00 , 8.21174e+00 , 8.62378e+00 , & 2.70834e-01 , 4.59350e-01 , 9.80330e-01 , 1.54856e+00 , & 2.14570e+00 , 2.75617e+00 , 3.36813e+00 , 3.97376e+00 , & 4.56847e+00 , 5.14997e+00 , 5.71743e+00 , 6.27093e+00 , & 6.81101e+00 , 7.33850e+00 , 7.85432e+00 , 8.35942e+00 , & 8.85471e+00 , 9.34106e+00 , 9.81925e+00 , 1.02900e+01 , & 3.99120e-01 , 6.46435e-01 , 1.31507e+00 , 2.03136e+00 , & 2.76957e+00 , 3.50946e+00 , 4.23832e+00 , 4.94958e+00 , & 5.64069e+00 , 6.31137e+00 , 6.96253e+00 , 7.59564e+00 , & 8.21234e+00 , 8.81426e+00 , 9.40295e+00 , 9.97983e+00 , & 1.05461e+01 , 1.11030e+01 , 1.16515e+01 , 1.21924e+01 , & 5.64538e-01 , 8.82405e-01 , 1.72870e+00 , 2.61770e+00 , & 3.51364e+00 , 4.39362e+00 , 5.24695e+00 , 6.07052e+00 , & 6.86499e+00 , 7.63269e+00 , 8.37647e+00 , 9.09919e+00 , & 9.80351e+00 , 1.04918e+01 , 1.11660e+01 , 1.18280e+01 , & 1.24792e+01 , 1.31209e+01 , 1.37543e+01 , 1.43803e+01 , & 7.72889e-01 , 1.17552e+00 , 2.23335e+00 , 3.31894e+00 , & 4.38704e+00 , 5.41672e+00 , 6.40295e+00 , 7.34788e+00 , & 8.25618e+00 , 9.13289e+00 , 9.98272e+00 , 1.08097e+01 , & 1.16174e+01 , 1.24086e+01 , 1.31857e+01 , 1.39506e+01 , & 1.47049e+01 , 1.54501e+01 , 1.61872e+01 , 1.69173e+01 & /) ,(/NTHICK,NFTAB/)) END IF DO I=1,NK DO J=1,NTHICK IF (SIG(I)*TPIINV.LT.SIS1FTABLE(1)) THEN SIS2ALPHAS(J,I) = SIS1ALPHATABLE(J,1) ELSE IF (SIG(I)*TPIINV.GT. SIS1FTABLE(NFTAB)) THEN SIS2ALPHAS(J,I) = SIS1ALPHATABLE(J,NFTAB) ELSE IND = 1 DO K = 1, NFTAB-1 IF (SIS1FTABLE(K).LT.SIG(I)*TPIINV) IND = K END DO X=(SIG(I)*TPIINV-SIS1FTABLE(IND))/(SIS1FTABLE(IND+1)-SIS1FTABLE(IND)) SIS2ALPHAS(J,I)=SIS1ALPHATABLE(J,IND)*(1-X)+SIS1ALPHATABLE(J,IND+1)*X END IF ! WRITE(998,*) I, J, SIG(I)*TPIINV,SIS1FTABLE(NFTAB), X, IND, SIS2ALPHAS(J,I),SIS1ALPHATABLE(J,NFTAB) END DO END DO ! ! SIS1ALPHATABLE2 = reshape((/ & 0.000001693306, 0.000001694535, 0.000001709985, 0.000001718371, 0.000001718689, & 0.000001715291, 0.000001711765, 0.000001709562, 0.000001708573, 0.000001708047, & 0.000001707216, 0.000001705508, 0.000001702489, 0.000001697757, 0.000001690926, & 0.000001681715, 0.000001670039, 0.000001655970, 0.000001639540, 0.000001620666, & 0.000002345944, 0.000002373942, 0.000002352788, 0.000002336058, 0.000002340157, & 0.000002361578, 0.000002389818, 0.000002414910, 0.000002430787, 0.000002435510, & 0.000002429671, 0.000002414451, 0.000002390516, 0.000002358094, 0.000002317788, & 0.000002271221, 0.000002220674, 0.000002167641, 0.000002111561, 0.000002051905, & 0.000003223057, 0.000003267049, 0.000003224977, 0.000003204685, 0.000003222416, & 0.000003261346, 0.000003299552, 0.000003321243, 0.000003319410, 0.000003293875, & 0.000003247777, 0.000003184750, 0.000003107742, 0.000003019246, 0.000002922032, & 0.000002819435, 0.000002714739, 0.000002609966, 0.000002505271, 0.000002400993, & 0.000004390462, 0.000004433012, 0.000004403102, 0.000004405652, 0.000004433604, & 0.000004455373, 0.000004447247, 0.000004400212, 0.000004316137, 0.000004201667, & 0.000004064075, 0.000003909816, 0.000003744650, 0.000003573915, 0.000003402281, & 0.000003233187, 0.000003068613, 0.000002909686, 0.000002757781, 0.000002614371, & 0.000005932907, 0.000005957402, 0.000005982223, 0.000006017129, 0.000006011597, & 0.000005931834, 0.000005775734, 0.000005561133, 0.000005309860, 0.000005038831, & 0.000004758982, 0.000004477929, 0.000004202314, 0.000003938073, 0.000003689194, & 0.000003456677, 0.000003239230, 0.000003035853, 0.000002848219, 0.000002677864, & 0.000007962387, 0.000007965182, 0.000008068497, 0.000008081519, 0.000007918056, & 0.000007586175, 0.000007143591, 0.000006652694, 0.000006157274, 0.000005680186, & 0.000005230854, 0.000004813199, 0.000004429817, 0.000004082375, 0.000003770294, & 0.000003490093, 0.000003236746, 0.000003006720, 0.000002800108, 0.000002616762, & 0.000010630745, 0.000010631503, 0.000010754695, 0.000010555056, 0.000009993969, & 0.000009206201, 0.000008341755, 0.000007500204, 0.000006728370, 0.000006039825, & 0.000005432987, 0.000004901489, 0.000004438118, 0.000004035401, 0.000003685333, & 0.000003379621, 0.000003110711, 0.000002873062, 0.000002663457, 0.000002478682, & 0.000014147503, 0.000014179502, 0.000014068022, 0.000013255509, 0.000011955952, & 0.000010517681, 0.000009161184, 0.000007970672, 0.000006956434, 0.000006102184, & 0.000005386202, 0.000004787198, 0.000004285040, 0.000003861250, 0.000003500220, & 0.000003190247, 0.000002923212, 0.000002692594, 0.000002491216, 0.000002312444, & 0.000018803447, 0.000018843681, 0.000017888469, 0.000015846997, 0.000013469781, & 0.000011291688, 0.000009484334, 0.000008030183, 0.000006862685, 0.000005921675, & 0.000005161810, 0.000004546889, 0.000004045234, 0.000003629478, 0.000003278778, & 0.000002980118, 0.000002726294, 0.000002510734, 0.000002323048, 0.000002154027, & 0.000024999293, 0.000024767047, 0.000021866213, 0.000017907275, 0.000014286707, & 0.000011453674, 0.000009339677, 0.000007754652, 0.000006540404, 0.000005593636, & 0.000004848620, 0.000004258144, 0.000003783667, 0.000003393700, 0.000003065741, & 0.000002787091, 0.000002551635, 0.000002353114, 0.000002179924, 0.000002022176, & 0.000033276027, 0.000031812088, 0.000025413911, 0.000019081783, 0.000014360014, & 0.000011113013, 0.000008878307, 0.000007287536, 0.000006111165, 0.000005216269, & 0.000004523259, 0.000003978625, 0.000003542278, 0.000003183888, 0.000002882748, & 0.000002627060, 0.000002410531, 0.000002226889, 0.000002066156, 0.000001920400, & 0.000044336466, 0.000039340417, 0.000027862847, 0.000019248910, 0.000013854771, & 0.000010491793, 0.000008293280, 0.000006775936, 0.000005680752, 0.000004862822, & 0.000004234499, 0.000003739067, 0.000003338786, 0.000003008719, 0.000002732485, & 0.000002498817, 0.000002299071, 0.000002126144, 0.000001974716, 0.000001841606, & 0.000059030780, 0.000046157443, 0.000028771759, 0.000018580063, 0.000013047269, & 0.000009812168, 0.000007739586, 0.000006322672, 0.000005313968, 0.000004570787, & 0.000004001034, 0.000003546046, 0.000003172357, 0.000002862735, 0.000002606289, & 0.000002391441, 0.000002205202, 0.000002038837, 0.000001893870, 0.000001773775, & 0.000078242681, 0.000050872198, 0.000028192943, 0.000017445281, 0.000012194039, & 0.000009215299, 0.000007293336, 0.000005968851, 0.000005030058, 0.000004344187, & 0.000003817214, 0.000003389423, 0.000003031921, 0.000002735020, 0.000002492733, & 0.000002292508, 0.000002116550, 0.000001954358, 0.000001814339, 0.000001706341, & 0.000102553944, 0.000052628108, 0.000026654291, 0.000016235235, 0.000011448810, & 0.000008738110, 0.000006951009, 0.000005700488, 0.000004811678, 0.000004163305, & 0.000003662892, 0.000003251855, 0.000002905015, 0.000002617308, 0.000002385271, & 0.000002195475, 0.000002027393, 0.000001869747, 0.000001734845, 0.000001635748, & 0.000131528803, 0.000051648919, 0.000024847696, 0.000015214276, 0.000010847813, & 0.000008336757, 0.000006656766, 0.000005471182, 0.000004622914, 0.000003999733, & 0.000003516266, 0.000003119051, 0.000002785020, 0.000002508708, 0.000002285704, & 0.000002102851, 0.000001941273, 0.000001790864, 0.000001662390, 0.000001566768, & 0.000162640998, 0.000049040647, 0.000023269162, 0.000014460197, 0.000010345625, & 0.000007943173, 0.000006352063, 0.000005238097, 0.000004434583, 0.000003834908, & 0.000003367217, 0.000002987629, 0.000002673539, 0.000002414062, 0.000002200682, & 0.000002022244, 0.000001866432, 0.000001726382, 0.000001605605, 0.000001507441, & 0.000190593962, 0.000046018523, 0.000022036249, 0.000013884361, 0.000009879926, & 0.000007534440, 0.000006029990, 0.000004998037, 0.000004246561, 0.000003673278, & 0.000003223583, 0.000002865466, 0.000002575963, 0.000002336404, 0.000002132787, & 0.000001957123, 0.000001806342, 0.000001677844, 0.000001564855, 0.000001460377, & 0.000208777708, 0.000043254960, 0.000020953050, 0.000013327645, 0.000009433705, & 0.000007161759, 0.000005746018, 0.000004790209, 0.000004085598, 0.000003536085, & 0.000003102390, 0.000002762488, 0.000002493229, 0.000002269788, 0.000002073862, & 0.000001900259, 0.000001753553, 0.000001634519, 0.000001527844, 0.000001418102, & 0.000213643676, 0.000040801760, 0.000019820045, 0.000012719494, 0.000009042103, & 0.000006888517, 0.000005551025, 0.000004645394, 0.000003968666, 0.000003433279, & 0.000003008733, 0.000002677964, 0.000002418184, 0.000002202481, 0.000002011120, & 0.000001840100, 0.000001696853, 0.000001583290, 0.000001480832, 0.000001370962, & 0.000208202335, 0.000038538552, 0.000018774435, 0.000012168047, 0.000008723021, & 0.000006681551, 0.000005399851, 0.000004523630, 0.000003863808, 0.000003339344, & 0.000002922742, 0.000002598268, 0.000002343776, 0.000002132781, 0.000001945837, & 0.000001779030, 0.000001639702, 0.000001529688, 0.000001430567, 0.000001324016, & 0.000198860865, 0.000036596393, 0.000018094381, 0.000011782298, 0.000008423914, & 0.000006432368, 0.000005194822, 0.000004356699, 0.000003726655, 0.000003224518, & 0.000002825485, 0.000002515982, 0.000002274472, 0.000002073961, 0.000001894595, & 0.000001733161, 0.000001598541, 0.000001493403, 0.000001398055, 0.000001292778, & 0.000188086831, 0.000035005432, 0.000017449034, 0.000011378985, 0.000008128438, & 0.000006206255, 0.000005018705, 0.000004216819, 0.000003612753, 0.000003129463, & 0.000002744757, 0.000002446723, 0.000002214548, 0.000002021421, 0.000001847651, & 0.000001690537, 0.000001559668, 0.000001458051, 0.000001365539, 0.000001262003, & 0.000176803511, 0.000033479878, 0.000016713766, 0.000010962586, 0.000007891418, & 0.000006057305, 0.000004907363, 0.000004122061, 0.000003528516, 0.000003054204, & 0.000002676573, 0.000002382977, 0.000002153342, 0.000001962633, 0.000001792438, & 0.000001639642, 0.000001512208, 0.000001412404, 0.000001322051, 0.000001223019, & 0.000157501278, 0.000031061181, 0.000013781133, 0.000009601049, 0.000007947466, & 0.000006760411, 0.000005655833, 0.000004680347, 0.000003892156, 0.000003277589, & 0.000002790517, 0.000002392986, 0.000002067558, 0.000001808761, 0.000001609460, & 0.000001454171, 0.000001323348, 0.000001206327, 0.000001111024, 0.000001045460, & 0.000007782240, 0.000002755736, 0.000000005024, 0.000000101927, 0.000088216814, & 0.009792270789, 0.030135654570, 0.009293514348, 0.001549599316, 0.000298014001, & 0.000055965309, 0.000007063938, 0.000000637549, 0.000000072256, 0.000000017854, & 0.000000008901, 0.000000004238, 0.000000001169, 0.000000000492, 0.000000001113 /), & (/NTHICK,NFTAB/)) DO I=1,NK DO J=1,NTHICK IF (SIG(I)*TPIINV.LT. SIS1FTABLE(1)) THEN SIS2ALPHA2(J,I) = SIS1ALPHATABLE2(J,1) ELSE IF (SIG(I)*TPIINV.GT. SIS1FTABLE(NFTAB)) THEN SIS2ALPHA2(J,I) = SIS1ALPHATABLE2(J,NFTAB) ELSE IND = 1 DO K = 1, NFTAB-1 IF (SIS1FTABLE(K).LT.SIG(I)*TPIINV) IND = K END DO X=(SIG(I)*TPIINV-SIS1FTABLE(IND))/(SIS1FTABLE(IND+1)-SIS1FTABLE(IND)) SIS2ALPHA2(J,I)=SIS1ALPHATABLE2(J,IND)*(1-X)+SIS1ALPHATABLE2(J,IND+1)*X END IF END DO END DO ! ! -------------------------------------------------------------------- / ! 2. Fills array of ICEDMAX to ICEDAVE ! DO I=1,NICED ICEDAVETAB(I) = W3FSD_DAVE(IS2PARS(9),REAL(I),IS2PARS(8)) ENDDO ! ! -------------------------------------------------------------------- / ! 2. Defines and diagonalizes the scattering matrix ! ALLOCATE(IS2SCATMAT(NTH,NTH)) ALLOCATE(IS2EIGVEC(NTH,NTH)) ALLOCATE(IS2EIGVAL(NTH)) ! DO I=1,NTH DO J=1,NTH ! This is for isotropic back-scatter ! IS2SCATMAT(I,J)=-1./DBLE(NTH) ! Other example that looks like figure 12 in Masson & LeBlond IS2SCATMAT(I,J)=-1./DBLE(NTH)*(2.*EC2(ABS(I-J)+1)**2+0.8*ECOS(ABS(I-J)+1)**3) IF (ECOS(ABS(I-J)+1).LT.0.001) IS2SCATMAT(I,J)=IS2SCATMAT(I,J)-1./DBLE(NTH)*0.8*ES2(ABS(I-J)+1) END DO !WRITE(997,'(36G16.8)') IS2SCATMAT(I,:) ! Now removes sum from diagonal to enforce energy conservation ... IS2SCATMAT(I,I)=IS2SCATMAT(I,I)-SUM(IS2SCATMAT(I,1:NTH)) END DO CALL DIAGONALIZE(IS2SCATMAT,IS2EIGVAL,IS2EIGVEC,nrot) DO I=1,NTH IS2EIGVAL(I)=MAX(0.d0,IS2EIGVAL(I)) !WRITE(994,'(36G16.8)') I,IS2EIGVAL(I) !WRITE(995,'(36G16.8)') IS2EIGVEC(I,:) !WRITE(996,'(36G16.8)') IS2SCATMAT(I,:) END DO !CLOSE(994) !CLOSE(995) !CLOSE(996) !CLOSE(997) RETURN !/ !/ End of INSIS2 ----------------------------------------------------- / !/ END SUBROUTINE INSIS2 !/ ------------------------------------------------------------------- / !> !> @brief Wave scattering in the MIZ, adapted from Dumont et al. !> !> @details This scattering routine allows the estimation of the !> maximum floe size and an estimate of the creep-induced dissipation. !> For the scattering, it is based on the normal incidence results of !> Kohout and Meylan which are provided in a table. !> !> @param[in] A Action density spectrum (1-D). !> @param[in] DEPTH Water depth. !> @param[in] CICE Sea ice concentration. !> @param[in] ICEH Ice thickness. !> @param[inout] ICEF Maximum floe size (updated). !> @param[in] ICEDMAX Maximum floe size. !> @param[in] IX Not used. !> @param[in] IY Not used. !> @param[out] S Source term (1-D version). !> @param[out] D Diagonal part of scattering (1-D version). !> @param[out] DISSIP Diagonal dissipation term (1-D version) !> @param[in] WN Wave number. !> @param[in] CG Group speed. !> @param[in] WN_R Wave number in ice. !> @param[in] CG_ICE Group speed in ice. !> @param[out] R Ratio of energy to wave energy without ice. !> !> @author P. Nicot !> @author F. Ardhuin !> @author G. Boutin !> @date 04-May-2016 !> SUBROUTINE W3SIS2 (A, DEPTH, CICE, ICEH, ICEF, ICEDMAX, IX, IY, & S, D, DISSIP, WN, CG, WN_R, CG_ICE, R) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | P. Nicot & F. Ardhuin & G. Boutin | !/ | FORTRAN 90 | !/ | Last update : 04-May-2016 | !/ +-----------------------------------+ !/ !/ 16-Mar-2014 : Origination. ( version 4.18 ) !/ 19-Sep-2014 : Correcting group speed factor ( version 5.03 ) !/ 20-Sep-2014 : Adding back-scattered energy ( version 5.03 ) !/ 27-Aug-2015 : Add breaking criterion, WIM1d ( version 5.05 ) !/ (ref. Williams, 2012) !/ 02-Nov-2015 : Integration of strain over bandwidth( version 5.05 ) !/ 13-Jan-2016 : Changed initialization of ICEDMAX ( version 5.10 ) !/ 06-Feb-2016 : Added IICEHMIN and creep dissipation( version 5.10 ) !/ 10-Mar-2016 : Added depth and call to Liu disp. ( version 5.10 ) !/ 02-May-2016 : Call to Liu disp moved to w3srce ( version 5.10 ) ! ! 1. Purpose : ! Wave scattering in the MIZ, adapted from Dumont et al. ! !/ ------------------------------------------------------------------- / ! ! 2. Method : ! This scattering routine allows the estimation of the maximum floe ! size and an estimate of the creep-induced dissipation. ! For the scattering, it is based on the normal incidence results of ! Kohout and Meylan which are provided in a table. ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! A R.A. I Action density spectrum (1-D) ! DEPTH Real I Water depth ! CICE Real I Sea ice concentration ! ICEH Real I ice thickness ! ICEF Real I/O Maximum floe size (updated) ! ICEDMAX Real I Maximum floe size ! 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) ! 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 ! CG_ICE R.A. I Group speed in ice ! R R.A. O Ratio of energy to wave energy without ice ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! W3SRCE Subr. W3SRCEMD Source term integration. ! W3EXPO Subr. N/A ASCII Point output post-processor. ! W3EXNC Subr. N/A NetCDF Point output post-processor. ! GXEXPO Subr. N/A GrADS point output post-processor. ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! If ice concentration is zero, no calculations are made. ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! !/S Enable subroutine tracing. ! !/T Enable general test output. ! 2-D print plot of source term. ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE W3ODATMD, ONLY: NDSE USE CONSTANTS, ONLY: TPIINV, PI, TPI, GRAV, DWAT USE W3SERVMD, ONLY: EXTCDE USE W3GDATMD, ONLY: NK, NTH, NSPEC, SIG, SIG2, DDEN, IS2PARS, XFR, & IICEHMIN,IICESMOOTH #ifdef W3_T USE W3ODATMD, ONLY: NDST #endif #ifdef W3_S USE W3SERVMD, ONLY: STRACE #endif #ifdef W3_T USE W3ARRYMD, ONLY: PRT2DS #endif USE W3DISPMD ! IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list REAL, INTENT(IN) :: A(NSPEC), DEPTH, CICE, ICEH, ICEDMAX INTEGER, INTENT(IN) :: IX, IY 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) !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ #ifdef W3_S INTEGER, SAVE :: IENT = 0 #endif INTEGER :: IK, IKP1, IKM1, ITH, ITH2, IS, IS2, IND1, IND2 REAL :: W INTEGER :: IKBREAK, NSUM LOGICAL :: BRK_CRIT_W(NK), BRK_CRIT REAL :: ALPHA, STRAIN_C, WAMP(NK), D_FLEX_FAIL, & SMOOTHD, TAU_D, S_D(NK), ALPHA_D, DELTA_D,B_COLE, & DMAX, S_ATT, FACTOR, BETA REAL :: ICEDAVE(NK), CURVTOSTRAIN, CREEPFAC, MP2, B, ICEF_CREEP REAL :: SUMALLDIR, SUMA, SUME, CURVSPEC(NK), ESPEC(NK),STRAIN REAL, PARAMETER :: YOUNG = 5.49E+9 ! Young modulus REAL, PARAMETER :: POISSON = 0.3 ! Poisson Ratio REAL :: SIGMA_C REAL, PARAMETER :: DENS = 1025.0 ! ice density REAL :: GAMMA_TOY REAL, DIMENSION(NK) :: WN_I, WN_RP, WSQ, WLG, WLG_I, CG_I, & CURV, CGRATIO, CG_EFF, DUMMY, ALPHA_DISP #ifdef W3_T REAL :: SOUT(NK,NTH) #endif !/ !/ ------------------------------------------------------------------- / !/ #ifdef W3_S CALL STRACE (IENT, 'W3SIS1') #endif ! ! 0. Initializations ------------------------------------------------ * ! S = 0. D = 0. DISSIP = 0. DUMMY = WN WN_I = 0. WN_RP = 0. CG_I = 0. CG_EFF = 0. SIGMA_C = IS2PARS(19) MP2=(1-POISSON**2) GAMMA_TOY = 2 + log(0.9)/log(2.) ! Variables from Cole et al. 1995 ALPHA_D=0.54 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(1:NK)*TAU_D) DELTA_D=IS2PARS(21) IF (IS2PARS(9).GT.0) ICEDMIN = IS2PARS(9) IF (IS2PARS(12).GT.0) THEN B=IS2PARS(12) ! 2 is the ratio of to (Hrms/2)^4 for a Rayleigh distribution ! 0.375 is the average of cos^4 ! 0.4 is 2/5 CREEPFAC = -2*(0.25/(IS2PARS(15)+2))*0.375*B*ICEH**(IS2PARS(15)+2) & *(YOUNG/(2*B*MP2))**(IS2PARS(15)+1)/(DWAT*GRAV) ELSE CREEPFAC=1. ENDIF WLG = TPI/WN ICEF_CREEP=ICEF DMAX = ICEF BRK_CRIT = .FALSE. NSUM = NINT(0.3/(XFR-1.)) ! STRAIN_C = SIGMA_C*MP2/YOUNG ! Minimum floe size that can break D_FLEX_FAIL = 0.5* ( (PI**4*YOUNG*ICEH**3)/( 48*DENS*GRAV*MP2 ) )**.25 ! Estimates mean floe diameter from max floe diameter IF (DMAX.GT.NICED) THEN ICEDAVE=DMAX ELSE IF (IS2PARS(9).GT.0) THEN ICEDAVE= ICEDAVETAB(MAX(1,NINT(DMAX))) ELSE DO IK=1,NK IF (IS2PARS(14)*(TPI/WN_R(IK)).LT.DMAX) THEN ICEDAVE(IK)=W3FSD_DAVE( (TPI/WN_R(IK))*IS2PARS(14),DMAX,IS2PARS(8)) ELSE ICEDAVE(IK)=DMAX ENDIF END DO ENDIF ENDIF ! #ifdef W3_T SOUT = 0. #endif ! IF (CICE .GT. 0) THEN ! ! 1. Calculate wavelength, Robinson-Palmer dispersion relation ! (should be tabulated ) ! Warning: it only applies to unbroken ice. ! IF (IS2PARS(6).GT.0.5) THEN CALL W3RPWNICE(ICEH,WN_I,WN_RP,CG_I) ELSE WN_I = WN_R CG_I = CG_ICE WN_RP=WN*0. END IF ! WLG_I = TPI/WN_I ! Ice wavelength WSQ = (WLG/WLG_I) ! IF (IS2PARS(16).GT.0.5) THEN CGRATIO(:)=CG_I(:)/CG(:) ELSE CGRATIO(:)=1 END IF ! ! 2. gets reflection coefficient from table ! IND1 = 1+FLOOR((ICEH-THICK1)/DTHICK) IND2 = IND1+1 ! defines weight for interpolation of ice thickness W = (ICEH-THICK1)/DTHICK - (IND1-1) IF (IND1.LT.1) THEN IND1 = 1 W = 0 ELSE IF (IND2.GT.NTHICK) THEN IND2 = NTHICK IND1 = NTHICK W = 0 END IF ! DO IK = 1,NK ! ! Spatial decay scale taken from table. This corresponds the the values shown ! in Dumont et al. (JGR 2011, fig 3). ! CG_EFF(IK) = CICE*CG_I(IK) + (1-CICE)*CG(IK) ! Note by FA: dissipation should be done by ICx not by ISx ! the RP damping is thus defined by an optional IS2PARS(7), which is 0 by default ALPHA = -1.*( (SIS2ALPHAS(IND1,IK)*(1-W)+SIS2ALPHAS(IND2,IK)*W)/ICEDAVE(IK) & +2.*IS2PARS(7)*WN_RP(IK) )*CICE ! ! Additional scattering for pack ice defined by IS2PARS(4:5) (see Squire et al. GRL 2009) ! ALPHA = IS2PARS(1) * ALPHA -2.*IS2PARS(4)*EXP(-1.*IS2PARS(5)/SIG(IK)) IF (IS2PARS(11).GT.0) THEN IF (CICE.LT.0.2) ALPHA = 0. IF (CICE.GT.0.8) ALPHA = 0. IF (CICE.GE.0.2.AND.CICE.LE.0.8) ALPHA = ALPHA*(CICE-0.2)*(0.8-CICE) END IF ! ! time decay ! BETA = ALPHA * CICE * CG_EFF(IK) ! ! 3. attenuation due to scattering for all spectral components ! with added backscattering for energy conservation ( if IS2PARS(2).EQ.1) ! SUMALLDIR= 0. SUMA = 0. CURVTOSTRAIN = (0.25*MAX(ICEH,IICEHMIN)**2) DO ITH = 1,NTH IS = ITH+(IK-1)*NTH D(IS) = BETA S(IS) = BETA * A(IS) SUMALLDIR = SUMALLDIR + S(IS) SUMA = SUMA + A(IS) END DO ! loop over directions ! ! R is the ratio of energy (including bending of ice) to wave energy without ice ! Wadhams 1973 eq. 34, warning, his ice thickness is 2*h ! Warning : R uses DMAX=ICEF, even if IS2DUPDATE=F ! IF (IICESMOOTH) THEN IF (IS2PARS(14)*WLG_I(IK).LT.DMAX) THEN SMOOTHD=TANH((DMAX-IS2PARS(14)*WLG_I(IK))/(DMAX*IS2PARS(13))) ELSE SMOOTHD=0. END IF ELSE SMOOTHD=1. END IF ! R(IK) =1+IS2PARS(16)*SMOOTHD*4*YOUNG*ICEH**3*(PI/WLG_I(IK))**4/(3*DWAT*GRAV*MP2) ! ! Converting action to surface elevation variance SUME with units m^2 ! SUME = SUMA*DDEN(IK) / CG(IK) / (R(IK)*CGRATIO(IK)) ! ! CURVSPEC is the curvature variance = elevation variance * k^4 ! CURVSPEC (IK) = SUME * (2*PI/ WLG_I(IK))**4 ESPEC (IK) = SUME SUMALLDIR = SUMALLDIR / REAL(NTH) ! ! Adds the scattered energy isotropically to conserve the energy ! This may not be a very good scheme numerically. Another possible ! approach is the matrix inversion used for bottom scattering (w3sbs1md.ftn) ! S(1+(IK-1)*NTH:IK*NTH)=S(1+(IK-1)*NTH:IK*NTH)-SUMALLDIR*IS2PARS(2) END DO ! loop over wavenumbers IK ! ! 4. update of floe size ! IF (IS2PARS(10).LT.0.5) THEN ! resets max floe size to the last forcing or initial value DMAX = ICEDMAX ICEF = ICEDMAX END IF ! DO IK = 1, NK ! CURV is the variance of the curvature integrated over a finite bandwidth CURV(IK) = SUM(CURVSPEC(MAX(1,IK-NSUM):MIN(NK,IK+NSUM))) ! Now converts curvature variance to strain variance END DO ! end of loop on IK ! ! If IS2PARS(3)=IS2BREAK is set to true in ww3_grid, then activates ice break-up ! IF (IS2PARS(3).GT.0.5) THEN IKBREAK=0 DO IK = 1, NK STRAIN = CURV(IK)*CURVTOSTRAIN IF (D_FLEX_FAIL .LT. DMAX) THEN ! Note that Williams et al. used IS2PARS(17)=SQRT(2), Here our default is 3.6 WAMP(IK)= IS2PARS(17)*SQRT(STRAIN) ! IF (IS2PARS(9).EQ.0) THEN ICEDMIN=(TPI/WN_R(IK))*IS2PARS(14) END IF BRK_CRIT_W(IK) = WAMP(IK) .GT. STRAIN_C .AND. WLG_I(IK)/2 .GT. ICEDMIN .AND. WLG_I(IK)/2 .LT. DMAX & .AND. WLG_I(IK)/2 .GT. D_FLEX_FAIL ! IF (BRK_CRIT_W(IK)) THEN IKBREAK=IK BRK_CRIT = .TRUE. END IF END IF END DO ! end of loop on IK ! ! 4.b Correction for bias introduced by the finite bandwidth sum .. . ! IF (BRK_CRIT) THEN DO IK=MAX(IKBREAK-NSUM,1),IKBREAK,1 ! Modified by F.A. on Jan. 31, 2017: uses the maximum of CURVSPEC instead of CURV. ! this is better for very narrow spectra. IF (CURVSPEC(IK).GE.CURVSPEC(IKBREAK).AND.DMAX.GE.(WLG_I(IK)/2)) THEN IKBREAK = IK END IF END DO ! DMAX = WLG_I(IKBREAK)/2 ! ! Uses a weighting by CURVSPEC to have a continuous shift of DMAX .. ! IKP1=MIN(IKBREAK+1,NK) IKM1=MAX(IKBREAK-1,1) IF (BRK_CRIT_W(IKP1).AND.BRK_CRIT_W(IKM1)) THEN DMAX = (WLG_I(IKBREAK)*CURVSPEC(IKBREAK) & +WLG_I(IKP1)*CURVSPEC(IKP1)+WLG_I(IKM1)*CURVSPEC(IKM1)) & /(2.*(CURVSPEC(IKBREAK)+CURVSPEC(IKM1)+CURVSPEC(IKP1))) END IF ! ICEF = DMAX END IF END IF !end of test (IS2PARS(3).GT.0.5) ! ! 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. ! 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 ! 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 ! 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+DISSIP*CICE*A END DO ! end of loop on IK ENDIF ! end of test (IS2PARS(12).GT.0) ! ! 6. Case of no scattering nor dissipation ! ELSE DMAX = 0. ICEF = 0. END IF ! end of test (CICE .GT. 0 .AND. ICEDAVE .GT. 0) ! #ifdef W3_T DO IK = 1, NK DO ITH = 1, NTH IS = ITH+(IK-1)*NTH SOUT(IK,ITH) = S(IS) END DO END DO CALL PRT2DS (NDST, NK, NK, NTH, SOUT, SIG(1:NK), ' ', 1., & 0.0, 0.001, 'Diag Sir1', ' ', 'NONAME') #endif ! ! Formats 8000 FORMAT (' TEST W3SIS2 : ALPHA :',E10.3) ! !/ !/ End of W3SIS2 ----------------------------------------------------- / !/ END SUBROUTINE W3SIS2 !/ ------------------------------------------------------------------- / !> !> @brief NA. !> !> @param[in] ICEH !> @param[inout] WN_I !> @param[inout] DAMPING !> @param[inout] CG_I !> !> @author C. Sevigny !> @date 27-Aug-2015 SUBROUTINE W3RPWNICE(ICEH,WN_I,DAMPING,CG_I) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | C. Sevigny | !/ | FORTRAN 90 | !/ | Last update : 27-Aug-2015 | !/ +-----------------------------------+ !/ !/ 27-Aug-2015 : Origination. ( version 5.10 ) ! ! 1. Purpose : ! !/ ------------------------------------------------------------------- / ! ! 2. Method : ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS, ONLY: GRAV USE W3GDATMD, ONLY: NK, SIG IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list REAL, INTENT(IN) :: ICEH REAL, INTENT(INOUT) :: WN_I(:), DAMPING(:), CG_I(:) !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ DOUBLE COMPLEX :: WN_ROOT, GS1 INTEGER :: IK REAL :: FLEX_RIGID REAL, PARAMETER :: VISC_RP = 10 REAL, PARAMETER :: DENS = 1025.0 REAL, PARAMETER :: DENS_ICE = 922.5 REAL, PARAMETER :: POISSON = 0.3 REAL, PARAMETER :: YOUNG = 5.49E+9 FLEX_RIGID = YOUNG * ICEH**3 /(12 *(1-POISSON**2) ) ! Guess value for roots GS1 = CMPLX(SIG(1)**2/GRAV,0.) DO IK=1,NK CALL FINDROOTS_NR(GS1,0.,WN_ROOT,ICEH,SIG(IK)) WN_I(IK) = REAL(WN_ROOT) CALL FINDROOTS_NR(GS1,VISC_RP,WN_ROOT,ICEH,SIG(IK)) DAMPING(IK) = AIMAG(WN_ROOT) GS1 = WN_I(IK) CG_I(IK) = (5* FLEX_RIGID*WN_I(IK)**4 + DENS*GRAV - DENS_ICE*ICEH*SIG(IK)**2) & /(2*SIG(IK)*(DENS+DENS_ICE*WN_I(IK)*ICEH)) END DO !/ END SUBROUTINE W3RPWNICE !/ ------------------------------------------------------------------- / ! SUBROUTINE FINDROOTS_NR(FUNCD,X0,VISC_RP,WN_ROOT) !> !> @brief NA. !> !> @param[in] GUESS !> @param[in] VISC_RP !> @param[in] ICEH !> @param[in] SIGMA !> @param[inout] X !> !> @author C. Sevigny !> @date 27-Aug-2015 !> SUBROUTINE FINDROOTS_NR(GUESS,VISC_RP,X,ICEH,SIGMA) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | C. Sevigny | !/ | FORTRAN 90 | !/ | Last update : 27-Aug-2015 | !/ +-----------------------------------+ !/ !/ 27-Aug-2015 : Origination. ( version 5.10 ) ! ! 1. Purpose : ! !/ ------------------------------------------------------------------- / ! ! 2. Method : ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE W3ODATMD, ONLY: NDSE USE W3SERVMD, ONLY: EXTCDE !/ !/ ------------------------------------------------------------------- / IMPLICIT NONE !/ ------------------------------------------------------------------- / !/ Parameter list REAL, INTENT(IN) :: VISC_RP, ICEH, SIGMA double COMPLEX, INTENT(IN) :: GUESS double COMPLEX, intent(INOUT) :: X !/ ------------------------------------------------------------------- / !/ Local parameters !/ double complex :: FVAL, FDERIV, DX, X0 INTEGER, PARAMETER :: MAXIT = 300 REAL, PARAMETER :: TOL = 1E-9 INTEGER :: J LOGICAL :: UNFINISHED_ROOTS X0 = GUESS UNFINISHED_ROOTS = .TRUE. J = 1 DO WHILE (J .LT. MAXIT .AND. UNFINISHED_ROOTS) FVAL = FUNCD_FVAL(X0,VISC_RP,ICEH,SIGMA) FDERIV = FUNCD_FDERIV(X0,VISC_RP,ICEH,SIGMA) DX = FVAL/FDERIV X = X0-DX X0 = X IF (ABS(DX) .GT. TOL .OR. ABS(FVAL) .GT. 1) THEN J = J+1 ELSE UNFINISHED_ROOTS = .FALSE. END IF IF (J .GT. MAXIT) THEN WRITE (NDSE,1000) CALL EXTCDE ( 1 ) END IF END DO ! ! Formats 1000 FORMAT (/' *** ERROR FINDROOTS_NR *** '/ & ' ROOT NOT CONVERGED'/) END SUBROUTINE FINDROOTS_NR !/ ------------------------------------------------------------------- / !> !> @brief Computes the mean flow size from the max floe size. !> !> @param ICEDMIN !> @param ICEDMAX !> @param FRAGILITY !> @returns W3FSD_DAVE !> !> @author C. Sevigny !> @author F. Ardhuin !> @date 06-Nov-2015 !> FUNCTION W3FSD_DAVE(ICEDMIN,ICEDMAX,FRAGILITY) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | C. Sevigny & F. Ardhuin | !/ | FORTRAN 90 | !/ | Last update : 6-Nov-2015 | !/ +-----------------------------------+ !/ !/ 27-Aug-2015 : Origination. ( version 5.10 ) !/ 6-Nov-2015 : Uses a continuous DMAX->DAVE function(version 5.10 ) ! ! 1. Purpose : Computes the mean flow size from the max floe size ! !/ ------------------------------------------------------------------- / ! ! 2. Method : ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ICEDMIN Real I Minimum floe diameter ! FRAGILITY Real I Parameter that gives the power in the FSD power law ! ICEDMAX R.A. I/O Maximum floe diameter ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! None ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! INSIS2 Proc. W3SIS2MD Initialisation of parameters for IS2 ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / ! IMPLICIT NONE !/ !/ ------------------------------------------------------------------- / !/ Parameter list REAL, INTENT(IN) :: ICEDMIN, ICEDMAX, FRAGILITY REAL :: W3FSD_DAVE, W3FSD_DAVE2 !/ ------------------------------------------------------------------- / !/ Local parameters !/ ! ! analytic solution, if the FSD is given by P(x) = x^(-1-gam) for icedmin <= x <= icedmax and 0 elsewhere ! REAL :: GAM, MR, R, DENOM, ICEDMAXL REAL, PARAMETER :: xi = 2 ! ICEDMAXL=MAX(ICEDMIN*1.01,ICEDMAX) GAM = 2 + log(FRAGILITY)/log(2.) IF (ICEDMIN.EQ.0) THEN W3FSD_DAVE =ICEDMAXL ELSE W3FSD_DAVE = GAM/(GAM-1) * ( (ICEDMAXL**(-GAM+1)-ICEDMIN**(-GAM+1))& /(ICEDMAXL**(-GAM)-ICEDMIN**(-GAM)) ) END IF ! ! Other possibility: analytical solution to Toyota algorithm (F. Arduin) ! MR=log(ICEDMAXL/ICEDMIN)/log(2.) ! R=(1-FRAGILITY*xi**2)/(1-FRAGILITY*xi) ! DENOM = (1-(FRAGILITY*xi**2)**(MR+1.)) ! W3FSD_DAVE2 = MAX(ICEDMAXL*R* (1-(FRAGILITY*xi)**(MR+1.))/DENOM,ICEDMIN) ! WRITE(991,*) ICEDMAX, GAM, W3FSD_DAVE, W3FSD_DAVE2 END FUNCTION W3FSD_DAVE ! !/ ------------------------------------------------------------------- / !> !> @brief NA. !> !> @param WN_GUESS !> @param VISC_RP !> @param ICEH !> @param SIGMA !> @returns FVAL !> !> @author C. Sevigny !> @date 27-Aug-2015 !> FUNCTION FUNCD_FVAL(WN_GUESS, VISC_RP,ICEH,SIGMA) RESULT(FVAL) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | C. Sevigny | !/ | FORTRAN 90 | !/ | Last update : 27-Aug-2015 | !/ +-----------------------------------+ !/ !/ 27-Aug-2015 : Origination. ( version 5.10 ) ! ! 1. Purpose : ! !/ ------------------------------------------------------------------- / ! ! 2. Method : ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS, ONLY: GRAV IMPLICIT NONE REAL, INTENT(IN) :: VISC_RP, ICEH, SIGMA DOUBLE COMPLEX, INTENT(IN) :: WN_GUESS !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ REAL :: ALP, FLEX_RIGID, GAMMA, C5, DRAFT, DENS_GRAV double complex :: C1, FVAL REAL, PARAMETER :: DENS = 1025.0 REAL, PARAMETER :: DENS_ICE = 922.5 REAL, PARAMETER :: POISSON = 0.3 ! Effective modulus or strain modulus [Pa], YOUNG REAL, PARAMETER :: YOUNG = 5.49E+9 !/ !/ ------------------------------------------------------------------- / ! Length below which flexural felure cannot occur, flexural rigidity FLEX_RIGID = YOUNG * ICEH**3 /(12 *(1-POISSON**2) ) DRAFT = (DENS_ICE/DENS)*ICEH DENS_GRAV = DENS*GRAV ALP = SIGMA**2/GRAV ! Artificial viscosity (Robinson & Palmer, 1990) GAMMA = SIGMA * VISC_RP/DENS_GRAV C1 = CMPLX(1 - ALP*DRAFT, - 1.*GAMMA) C5 = FLEX_RIGID/DENS_GRAV FVAL = C5*WN_GUESS**5 + C1*WN_GUESS - ALP END FUNCTION FUNCD_FVAL !/ ------------------------------------------------------------------- / ! FUNCTION FUNCD(WN_GUESS, VISC_RP,ICEH) !> !> @brief NA. !> !> @param WN_GUESS !> @param VISC_RP !> @param ICEH !> @param SIGMA !> @returns FDERIV !> !> @author C. Sevigny !> @date 27-Aug-2015 !> FUNCTION FUNCD_FDERIV(WN_GUESS, VISC_RP,ICEH,SIGMA) RESULT(FDERIV) !/ !/ +-----------------------------------+ !/ | WAVEWATCH III NOAA/NCEP | !/ | C. Sevigny | !/ | FORTRAN 90 | !/ | Last update : 27-Aug-2015 | !/ +-----------------------------------+ !/ !/ 27-Aug-2015 : Origination. ( version 5.10 ) ! ! 1. Purpose : ! !/ ------------------------------------------------------------------- / ! ! 2. Method : ! ! 3. Parameters : ! ! Parameter list ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 4. Subroutines used : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 5. Called by : ! ! Name Type Module Description ! ---------------------------------------------------------------- ! ---------------------------------------------------------------- ! ! 6. Error messages : ! ! None. ! ! 7. Remarks : ! ! 8. Structure : ! ! See source code. ! ! 9. Switches : ! ! 10. Source code : ! !/ ------------------------------------------------------------------- / USE CONSTANTS, ONLY: GRAV IMPLICIT NONE REAL, INTENT(IN) :: VISC_RP, ICEH, SIGMA DOUBLE COMPLEX, INTENT(IN) :: WN_GUESS !/ !/ ------------------------------------------------------------------- / !/ Local parameters !/ REAL :: ALP, FLEX_RIGID, GAMMA, C5, DRAFT, DENS_GRAV double complex :: C1, FDERIV REAL, PARAMETER :: DENS = 1025.0 REAL, PARAMETER :: DENS_ICE = 922.5 REAL, PARAMETER :: POISSON = 0.3 REAL, PARAMETER :: YOUNG = 5.49E+9 !/ !/ ------------------------------------------------------------------- / ! ! Length below which flexural felure cannot occur, flexural rigidity ! FLEX_RIGID = YOUNG * ICEH**3 /(12 * (1-POISSON**2)) DRAFT = (DENS_ICE/DENS)*ICEH DENS_GRAV = DENS*GRAV ALP = SIGMA**2/GRAV ! ! Artificial viscosity (Robinson & Palmer, 1990) ! GAMMA = SIGMA * VISC_RP/DENS_GRAV C1 = CMPLX(1 - ALP*DRAFT, -1.*GAMMA) C5 = FLEX_RIGID/DENS_GRAV ! FDERIV = 5*C5*WN_GUESS**4 +C1 ! END FUNCTION FUNCD_FDERIV !/ !/ End of module W3SIS1MD -------------------------------------------- / !/ END MODULE W3SIS2MD