-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathmain.f95
125 lines (94 loc) · 2.05 KB
/
main.f95
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
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
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
PROGRAM wind
USE param
USE sub
! local parameters
REAL :: time
INTEGER :: n,ntot,nout
mode = 1
!**********
CALL INIT ! initialisation
!**********
! initialisation of Eulerian tracer
DO j = 0,ny+1
DO k = 0,nx+1
T(j,k) = 0.
TN(j,k) = 0.
END DO
END DO
! output of initial eta distribution
OPEN(10,file ='eta0.dat',form='formatted')
DO j = 1,ny
WRITE(10,'(253F12.6)')(eta(j,k),k=1,nx)
END DO
CLOSE(10)
! output of initial layer thickness distribution
OPEN(10,file ='h0.dat',form='formatted')
DO j = 1,ny
WRITE(10,'(253F12.6)')(hzero(j,k),k=1,nx)
END DO
CLOSE(10)
! set local parameters
! set epsilon for Shapiro filter
eps = 0.0
! runtime parameters
ntot = INT(2.*24.*3600./dt)
! output parameter
nout = INT(900./dt)
! open files for output
OPEN(10,file ='eta.dat',form='formatted')
OPEN(20,file ='h.dat',form='formatted')
OPEN(30,file ='u.dat',form='formatted')
OPEN(40,file ='v.dat',form='formatted')
OPEN(50,file ='T.dat',form='formatted')
!---------------------------
! simulation loop
!---------------------------
DO n = 1,ntot
time = REAL(n)*dt
write(6,*)"time (days)", time/(24.*3600.)
taux = 0.2!*MIN(time/(1.*24.*3600.),1.0)
tauy = 0.0
DO j = 1,ny
T(j,0) = 0.0
T(J,1) = 0.0
END DO
DO j = 48,60
T(j,0) = 10.0
T(J,1) = 10.0
END DO
! call predictor
CALL dyn
! updating including Shapiro filter
CALL shapiro
! cyclic boundary conditions
DO j = 1,ny
DO k = 1,nx
etan(j,0) = etan(j,nx)
etan(j,nx+1) = etan(j,1)
etan(0,k) = etan(ny,k)
etan(ny+1,k) = etan(1,k)
END DO
END DO
DO j = 0,ny+1
DO k = 0,nx+1
h(j,k) = hzero(j,k) + eta(j,k)
wet(j,k) = 1
IF(h(j,k)<hmin)wet(j,k) = 0
u(j,k) = un(j,k)
v(j,k) = vn(j,k)
T(j,k) = TN(j,k)
END DO
END DO
! data output
IF(MOD(n,nout)==0)THEN
DO j = 1,ny
WRITE(10,'(253F12.6)')(eta(j,k),k=1,nx)
WRITE(20,'(253F12.6)')(h(j,k),k=1,nx)
WRITE(30,'(253F12.6)')(u(j,k),k=1,nx)
WRITE(40,'(253F12.6)')(v(j,k),k=1,nx)
WRITE(50,'(253F12.6)')(T(j,k),k=1,nx)
END DO
WRITE(6,*)"Data output at time = ",time/(24.*3600.)
ENDIF
END DO ! end of iteration loop
END PROGRAM wind