Blame view

src/floe_breaking.f90 3.71 KB
81dede1c   Jérémy Baudry   first commit
1

4e5b4414   Jérémy Baudry   ajout commentaire...
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
!_______________________________________________________________________________

                        !DESCRIPTION: Here is the parametrization proposed by
                        !Williams et al. (2013)  for wave induced
                        !floe breaking. In this routine, the maximum floe size
                        !is computed using the significant strain imposed by
                        !passing waves in ice.


!_______________________________________________________________________________



                        !INTERFACE:
                        subroutine floe_breaking

                        !USES:
                        use parameters

                        !LOCAL PARAMETERS:
81dede1c   Jérémy Baudry   first commit
22

4e5b4414   Jérémy Baudry   ajout commentaire...
23
                        implicit none
4508a952   Jérémy Baudry   new version
24

4e5b4414   Jérémy Baudry   ajout commentaire...
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
                        double precision, allocatable :: WN_ice(:)
                        double precision, allocatable :: fnt(:)
                        double precision, allocatable :: k(:)
                        double precision              :: rho=1025
                        double precision              :: d
                        double precision              :: F
                        double precision              :: Y=5.5
                        double precision              :: v=0.3
                        integer                       :: acc=10000
                        double precision, allocatable :: W(:)
                        double precision, allocatable :: fnt_d(:)

                        double precision, allocatable :: S(:)
                        double precision              :: SS
                        double precision              :: WN_ice_d
                        double precision              :: RS=7.05e-05
                        double precision              :: m0_s
                        double precision              :: m0_a
                        double precision              :: m2_a
                        double precision              :: Tw
                        double precision              :: om_d
                        double precision              :: wl_w

                        allocate(fnt(acc+1))
                        allocate(k(acc+1))
                        allocate(WN_ice(nfreq))
                        allocate(W(nfreq))
                        allocate(fnt_d(acc+1))
                        allocate(S(nfreq))

!_______________________________________________________________________________

                F=Y*h(i)**3/12d0*(1-v**2)
                d=0.9*h(i)

                jj=1            !create a wavelength vector (this is used to
                                !find the root of the polynom fnt
                do ii=0,acc
                k(jj)=real(ii)*10d0/real(acc)
                jj=jj+1
                end do


                do j=1,nfreq
                fnt=F*k**5+rho*(g-d*omega(j)**2)*k-rho*omega(j)**2
                WN_ice(j)=k(minloc(abs(fnt),DIM=1))!calculate wavelength in ice
                
                end do

 
                W=g*WN_ice/omega**2
                S=(0.5*h(i))*WN_ice**2*W
                m0_s=sum(E(n,i,1:nfreq)*S**2)*domega !0th moment of the strain spectrum
                SS=2*sqrt(m0_s)                      !significant strain
                if (SS.gt.RS) then
81dede1c   Jérémy Baudry   first commit
80

81dede1c   Jérémy Baudry   first commit
81

4e5b4414   Jérémy Baudry   ajout commentaire...
82
83
84
85
                m0_a=sum(E(n,i,1:nfreq)*W**2)*domega 
                m2_a=sum(omega**2*E(n,i,1:nfreq)*W**2)*domega
                Tw=2*pi*sqrt(m0_a/m2_a)
                om_d=2*pi/Tw
81dede1c   Jérémy Baudry   first commit
86

4e5b4414   Jérémy Baudry   ajout commentaire...
87
88
                fnt_d=F*k**5+rho*(g-d*om_d**2)*k-rho*om_d**2
                WN_ice_d=k(minloc(abs(fnt_d),DIM=1))
4508a952   Jérémy Baudry   new version
89

4e5b4414   Jérémy Baudry   ajout commentaire...
90
91
                wl_w=2*pi/WN_ice_d
                Dmax(i)=max(Dmin,min(0.5*wl_w,Dmax(i)))
81dede1c   Jérémy Baudry   first commit
92

81dede1c   Jérémy Baudry   first commit
93

4e5b4414   Jérémy Baudry   ajout commentaire...
94
                end if
81dede1c   Jérémy Baudry   first commit
95
96
97



4e5b4414   Jérémy Baudry   ajout commentaire...
98
!_________________________________________________________________________________
81dede1c   Jérémy Baudry   first commit
99
end subroutine floe_breaking