source: codes/icosagcm/trunk/src/timeloop_gcm.f90 @ 12

Last change on this file since 12 was 12, checked in by ymipsl, 12 years ago

dynamico tree creation

YM
File size: 7.3 KB
Line 
1MODULE timeloop_gcm_mod
2  USE genmod
3  USE transfert_mod
4  USE etat0_mod
5 
6  INTEGER,PARAMETER :: euler=1, leapfrog=2, leapfrog_matsuno=3, adam_bashforth=4
7 
8CONTAINS
9 
10  SUBROUTINE timeloop
11  USE field_mod
12  USE domain_mod
13  USE write_field
14  USE dimensions
15  USE geometry
16  USE transfert_mod
17  USE metric
18  USE dissip_mod
19  USE ioipsl
20  USE caldyn_gcm_mod
21  USE theta2theta_rhodz_mod
22  USE etat0_mod
23  IMPLICIT NONE
24  TYPE(t_field),POINTER :: f_phis(:)
25  TYPE(t_field),POINTER :: f_theta(:)
26  TYPE(t_field),POINTER :: f_ps(:),f_psm1(:), f_psm2(:)
27  TYPE(t_field),POINTER :: f_u(:),f_um1(:),f_um2(:)
28  TYPE(t_field),POINTER :: f_theta_rhodz(:),f_theta_rhodzm1(:),f_theta_rhodzm2(:)
29  TYPE(t_field),POINTER :: f_dps(:),f_dpsm1(:), f_dpsm2(:)
30  TYPE(t_field),POINTER :: f_du(:),f_dum1(:),f_dum2(:)
31  TYPE(t_field),POINTER :: f_dtheta_rhodz(:),f_dtheta_rhodzm1(:),f_dtheta_rhodzm2(:)
32
33  REAL(rstd),POINTER :: phis(:)
34  REAL(rstd),POINTER :: ps(:) ,psm1(:), psm2(:)
35  REAL(rstd),POINTER :: u(:,:) , um1(:,:), um2(:,:)
36  REAL(rstd),POINTER :: theta_rhodz(:,:) , theta_rhodzm1(:,:), theta_rhodzm2(:,:)
37  REAL(rstd),POINTER :: dps(:), dpsm1(:), dpsm2(:)
38  REAL(rstd),POINTER :: du(:,:), dum1(:,:), dum2(:,:)
39  REAL(rstd),POINTER :: dtheta_rhodz(:,:),dtheta_rhodzm1(:,:),dtheta_rhodzm2(:,:)
40  REAL(rstd) :: dt
41  INTEGER :: ind
42  INTEGER :: it,i,j,n
43  INTEGER :: scheme
44  INTEGER :: matsuno_period
45  INTEGER :: itaumax
46 
47
48  dt=90.
49  CALL getin('dt',dt)
50
51  itaumax=100
52  CALL getin('itaumax',itaumax)
53 
54  scheme=leapfrog_matsuno
55  CALL getin('scheme',scheme)
56 
57  matsuno_period=5
58  CALL getin('matsuno_period',matsuno_period)
59  IF (scheme==leapfrog) matsuno_period=itaumax+1
60 
61
62
63  CALL allocate_field(f_phis,field_t,type_real)
64 
65  CALL allocate_field(f_ps,field_t,type_real)
66  CALL allocate_field(f_psm1,field_t,type_real)
67  CALL allocate_field(f_psm2,field_t,type_real)
68  CALL allocate_field(f_dps,field_t,type_real)
69  CALL allocate_field(f_dpsm1,field_t,type_real)
70  CALL allocate_field(f_dpsm2,field_t,type_real)
71
72  CALL allocate_field(f_u,field_u,type_real,llm)
73  CALL allocate_field(f_um1,field_u,type_real,llm)
74  CALL allocate_field(f_um2,field_u,type_real,llm)
75  CALL allocate_field(f_du,field_u,type_real,llm)
76  CALL allocate_field(f_dum1,field_u,type_real,llm)
77  CALL allocate_field(f_dum2,field_u,type_real,llm)
78
79  CALL allocate_field(f_theta_rhodz,field_t,type_real,llm)
80  CALL allocate_field(f_theta_rhodzm1,field_t,type_real,llm)
81  CALL allocate_field(f_theta_rhodzm2,field_t,type_real,llm)
82  CALL allocate_field(f_dtheta_rhodz,field_t,type_real,llm)
83  CALL allocate_field(f_dtheta_rhodzm1,field_t,type_real,llm)
84  CALL allocate_field(f_dtheta_rhodzm2,field_t,type_real,llm)
85
86  CALL allocate_caldyn
87 
88!  CALL etat0_academic(f_ps,f_phis,f_theta_rhodz,f_u)
89  CALL etat0_jablonowsky06(f_ps,f_phis,f_theta_rhodz,f_u)
90!  CALL test_etat0_jablonowsky06
91 
92  DO it=0,itaumax
93    PRINT *,"It No :",It
94
95    CALL caldyn(f_phis,f_ps,f_theta_rhodz,f_u, f_dps, f_dtheta_rhodz, f_du)
96
97    IF (scheme==Euler) THEN
98      CALL  euler_scheme
99    ELSE IF (scheme==leapfrog) THEN
100      CALL leapfrog_scheme
101    ELSE IF (scheme==leapfrog_matsuno) THEN
102      CALL  leapfrog_matsuno_scheme
103    ELSE IF (scheme==adam_bashforth) THEN
104      CALL adam_bashforth_scheme
105    ENDIF
106
107  ENDDO
108 
109  CONTAINS
110
111    SUBROUTINE Euler_scheme
112    IMPLICIT NONE
113      INTEGER :: ind
114     
115      DO ind=1,ndomain
116        ps=f_ps(ind) ; u=f_u(ind) ; theta_rhodz=f_theta_rhodz(ind)
117        dps=f_dps(ind) ; du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind)
118        ps(:)=ps(:)+dt*dps(:)
119        u(:,:)=u(:,:)+dt*du(:,:)
120        theta_rhodz(:,:)=theta_rhodz(:,:)+dt*dtheta_rhodz(:,:)
121      ENDDO
122     
123    END SUBROUTINE Euler_scheme
124   
125    SUBROUTINE leapfrog_scheme
126    IMPLICIT NONE
127      INTEGER :: ind
128
129      DO ind=1,ndomain
130        ps=f_ps(ind)   ; u=f_u(ind)   ; theta_rhodz=f_theta_rhodz(ind)
131        psm1=f_psm1(ind) ; um1=f_um1(ind) ; theta_rhodzm1=f_theta_rhodzm1(ind)
132        psm2=f_psm2(ind) ; um2=f_um2(ind) ; theta_rhodzm2=f_theta_rhodzm2(ind)
133        dps=f_dps(ind) ; du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind)
134
135        IF (it==0) THEN
136          psm2(:)=ps(:) ; theta_rhodzm2(:,:)=theta_rhodz(:,:) ; um2(:,:)=u(:,:) 
137
138          ps(:)=ps(:)+dt*dps(:)
139          u(:,:)=u(:,:)+dt*du(:,:)
140          theta_rhodz(:,:)=theta_rhodz(:,:)+dt*dtheta_rhodz(:,:)
141
142          psm1(:)=ps(:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:) ; um1(:,:)=u(:,:) 
143        ELSE
144       
145          ps(:)=psm2(:)+2*dt*dps(:)
146          u(:,:)=um2(:,:)+2*dt*du(:,:)
147          theta_rhodz(:,:)=theta_rhodzm2(:,:)+2*dt*dtheta_rhodz(:,:)
148
149          psm2(:)=psm1(:) ; theta_rhodzm2(:,:)=theta_rhodzm1(:,:) ; um2(:,:)=um1(:,:) 
150          psm1(:)=ps(:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:) ; um1(:,:)=u(:,:) 
151        ENDIF
152         
153      ENDDO
154    END SUBROUTINE leapfrog_scheme 
155 
156    SUBROUTINE leapfrog_matsuno_scheme
157    IMPLICIT NONE
158      INTEGER :: ind
159
160      DO ind=1,ndomain
161        ps=f_ps(ind)   ; u=f_u(ind)   ; theta_rhodz=f_theta_rhodz(ind)
162        psm1=f_psm1(ind) ; um1=f_um1(ind) ; theta_rhodzm1=f_theta_rhodzm1(ind)
163        psm2=f_psm2(ind) ; um2=f_um2(ind) ; theta_rhodzm2=f_theta_rhodzm2(ind)
164        dps=f_dps(ind) ; du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind)
165
166       
167        IF (MOD(it,matsuno_period+1)==0) THEN
168          psm1(:)=ps(:) ; um1(:,:)=u(:,:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:)
169          ps(:)=psm1(:)+dt*dps(:)
170          u(:,:)=um1(:,:)+dt*du(:,:)
171          theta_rhodz(:,:)=theta_rhodzm1(:,:)+dt*dtheta_rhodz(:,:)
172
173        ELSE IF (MOD(it,matsuno_period+1)==1) THEN
174
175          ps(:)=psm1(:)+dt*dps(:)
176          u(:,:)=um1(:,:)+dt*du(:,:)
177          theta_rhodz(:,:)=theta_rhodzm1(:,:)+dt*dtheta_rhodz(:,:)
178
179          psm2(:)=psm1(:) ; theta_rhodzm2(:,:)=theta_rhodzm1(:,:) ; um2(:,:)=um1(:,:) 
180          psm1(:)=ps(:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:) ; um1(:,:)=u(:,:) 
181
182        ELSE
183
184          ps(:)=psm2(:)+2*dt*dps(:)
185          u(:,:)=um2(:,:)+2*dt*du(:,:)
186          theta_rhodz(:,:)=theta_rhodzm2(:,:)+2*dt*dtheta_rhodz(:,:)
187
188          psm2(:)=psm1(:) ; theta_rhodzm2(:,:)=theta_rhodzm1(:,:) ; um2(:,:)=um1(:,:) 
189          psm1(:)=ps(:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:) ; um1(:,:)=u(:,:) 
190
191        ENDIF
192     
193      ENDDO
194     
195    END SUBROUTINE leapfrog_matsuno_scheme 
196         
197    SUBROUTINE adam_bashforth_scheme
198    IMPLICIT NONE
199      INTEGER :: ind
200
201      DO ind=1,ndomain
202        ps=f_ps(ind)   ; u=f_u(ind)   ; theta_rhodz=f_theta_rhodz(ind)
203        dps=f_dps(ind) ; du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind)
204        dpsm1=f_dpsm1(ind) ; dum1=f_dum1(ind) ; dtheta_rhodzm1=f_dtheta_rhodzm1(ind)
205        dpsm2=f_dpsm2(ind) ; dum2=f_dum2(ind) ; dtheta_rhodzm2=f_dtheta_rhodzm2(ind)
206
207        IF (it==0) THEN
208          dpsm1(:)=dps(:) ; dum1(:,:)=du(:,:) ; dtheta_rhodzm1(:,:)=dtheta_rhodz(:,:)
209          dpsm2(:)=dpsm1(:) ; dum2(:,:)=dum1(:,:) ; dtheta_rhodzm2(:,:)=dtheta_rhodzm1(:,:)
210        ENDIF
211             
212        ps(:)=ps(:)+dt*(23*dps(:)-16*dpsm1(:)+5*dpsm2(:))/12
213        u(:,:)=u(:,:)+dt*(23*du(:,:)-16*dum1(:,:)+5*dum2(:,:))/12
214        theta_rhodz(:,:)=theta_rhodz(:,:)+dt*(23*dtheta_rhodz(:,:)-16*dtheta_rhodzm1(:,:)+5*dtheta_rhodzm2(:,:))/12
215
216        dpsm2(:)=dpsm1(:) ; dum2(:,:)=dum1(:,:) ; dtheta_rhodzm2(:,:)=dtheta_rhodzm1(:,:)
217        dpsm1(:)=dps(:) ; dum1(:,:)=du(:,:) ; dtheta_rhodzm1(:,:)=dtheta_rhodz(:,:)
218
219      ENDDO
220     
221     
222    END SUBROUTINE adam_bashforth_scheme
223
224  END SUBROUTINE timeloop
225   
226   
227 
228END MODULE timeloop_gcm_mod
Note: See TracBrowser for help on using the repository browser.