Context Navigation

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

Last change on this file since 276 was 266, checked in by ymipsl, 10 years ago

Synchronize trunk and Saturn branch.
Merge modification from Saturn branch to trunk

YM

File size: 19.2 KB

Line
1	MODULE timeloop_gcm_mod
2	USE transfert_mod
3	USE icosa
4	PRIVATE
5
6	PUBLIC :: init_timeloop, timeloop
7
8	INTEGER, PARAMETER :: euler=1, rk4=2, mlf=3
9	INTEGER, PARAMETER :: itau_sync=10
10
11	TYPE(t_message),SAVE :: req_ps0, req_mass0, req_theta_rhodz0, req_u0, req_q0
12
13	TYPE(t_field),POINTER,SAVE :: f_q(:)
14	TYPE(t_field),POINTER,SAVE :: f_rhodz(:), f_mass(:), f_massm1(:), f_massm2(:), f_dmass(:)
15	TYPE(t_field),POINTER,SAVE :: f_phis(:), f_ps(:),f_psm1(:), f_psm2(:), f_dps(:)
16	TYPE(t_field),POINTER,SAVE :: f_u(:),f_um1(:),f_um2(:), f_du(:)
17	TYPE(t_field),POINTER,SAVE :: f_theta_rhodz(:),f_theta_rhodzm1(:),f_theta_rhodzm2(:), f_dtheta_rhodz(:)
18	TYPE(t_field),POINTER,SAVE :: f_hflux(:), f_wflux(:), f_hfluxt(:), f_wfluxt(:)
19
20	INTEGER,SAVE :: nb_stage, matsuno_period, scheme
21	!$OMP THREADPRIVATE(nb_stage, matsuno_period, scheme)
22
23	CONTAINS
24
25	SUBROUTINE init_timeloop
26	USE icosa
27	USE dissip_gcm_mod
28	USE caldyn_mod
29	USE etat0_mod
30	USE disvert_mod
31	USE guided_mod
32	USE advect_tracer_mod
33	USE physics_mod
34	USE mpipara
35	USE omp_para
36	USE trace
37	USE transfert_mod
38	USE check_conserve_mod
39	USE output_field_mod
40	USE write_field
41	IMPLICIT NONE
42
43	CHARACTER(len=255) :: def
44
45
46	IF (xios_output) itau_out=1
47	IF (.NOT. enable_io) itau_out=HUGE(itau_out)
48
49	! Time-independant orography
50	CALL allocate_field(f_phis,field_t,type_real,name='phis')
51	! Trends
52	CALL allocate_field(f_du,field_u,type_real,llm,name='du')
53	CALL allocate_field(f_dtheta_rhodz,field_t,type_real,llm,name='dtheta_rhodz')
54	! Model state at current time step (RK/MLF/Euler)
55	CALL allocate_field(f_ps,field_t,type_real, name='ps')
56	CALL allocate_field(f_mass,field_t,type_real,llm,name='mass')
57	CALL allocate_field(f_u,field_u,type_real,llm,name='u')
58	CALL allocate_field(f_theta_rhodz,field_t,type_real,llm,name='theta_rhodz')
59	! Model state at previous time step (RK/MLF)
60	CALL allocate_field(f_um1,field_u,type_real,llm,name='um1')
61	CALL allocate_field(f_theta_rhodzm1,field_t,type_real,llm,name='theta_rhodzm1')
62	! Tracers
63	CALL allocate_field(f_q,field_t,type_real,llm,nqtot,'q')
64	CALL allocate_field(f_rhodz,field_t,type_real,llm,name='rhodz')
65	! Mass fluxes
66	CALL allocate_field(f_hflux,field_u,type_real,llm) ! instantaneous mass fluxes
67	CALL allocate_field(f_hfluxt,field_u,type_real,llm) ! mass "fluxes" accumulated in time
68	CALL allocate_field(f_wflux,field_t,type_real,llm+1) ! vertical mass fluxes
69	CALL allocate_field(f_dmass,field_t,type_real,llm, name='dmass')
70
71	IF(caldyn_eta == eta_mass) THEN ! eta = mass coordinate (default)
72	CALL allocate_field(f_dps,field_t,type_real,name='dps')
73	CALL allocate_field(f_psm1,field_t,type_real,name='psm1')
74	CALL allocate_field(f_wfluxt,field_t,type_real,llm+1,name='wfluxt')
75	! the following are unused but must point to something
76	! f_massm1 => f_mass
77	ELSE ! eta = Lagrangian vertical coordinate
78	CALL allocate_field(f_massm1,field_t,type_real,llm, name='massm1')
79	! the following are unused but must point to something
80	f_wfluxt => f_wflux
81	f_dps => f_phis
82	f_psm1 => f_phis
83	END IF
84
85	def='runge_kutta'
86	CALL getin('scheme',def)
87
88	SELECT CASE (TRIM(def))
89	CASE('euler')
90	scheme=euler
91	nb_stage=1
92	CASE ('runge_kutta')
93	scheme=rk4
94	nb_stage=4
95	CASE ('leapfrog_matsuno')
96	scheme=mlf
97	matsuno_period=5
98	CALL getin('matsuno_period',matsuno_period)
99	nb_stage=matsuno_period+1
100	! Model state 2 time steps ago (MLF)
101	CALL allocate_field(f_theta_rhodzm2,field_t,type_real,llm)
102	CALL allocate_field(f_um2,field_u,type_real,llm)
103	IF(caldyn_eta == eta_mass) THEN ! eta = mass coordinate (default)
104	CALL allocate_field(f_psm2,field_t,type_real)
105	! the following are unused but must point to something
106	f_massm2 => f_mass
107	ELSE ! eta = Lagrangian vertical coordinate
108	CALL allocate_field(f_massm2,field_t,type_real,llm)
109	! the following are unused but must point to something
110	f_psm2 => f_phis
111	END IF
112
113	CASE default
114	PRINT*,'Bad selector for variable scheme : <', TRIM(def), &
115	' > options are <euler>, <runge_kutta>, <leapfrog_matsuno>'
116	STOP
117	END SELECT
118
119
120	CALL init_dissip
121	CALL init_caldyn
122	CALL init_guided
123	CALL init_advect_tracer
124	CALL init_check_conserve
125	CALL init_physics
126
127	CALL etat0(f_ps,f_mass,f_phis,f_theta_rhodz,f_u, f_q)
128
129	CALL transfert_request(f_phis,req_i0)
130	CALL transfert_request(f_phis,req_i1)
131	CALL writefield("phis",f_phis,once=.TRUE.)
132
133	CALL init_message(f_ps,req_i0,req_ps0)
134	CALL init_message(f_mass,req_i0,req_mass0)
135	CALL init_message(f_theta_rhodz,req_i0,req_theta_rhodz0)
136	CALL init_message(f_u,req_e0_vect,req_u0)
137	CALL init_message(f_q,req_i0,req_q0)
138
139	END SUBROUTINE init_timeloop
140
141	SUBROUTINE timeloop
142	USE icosa
143	USE dissip_gcm_mod
144	USE disvert_mod
145	USE caldyn_mod
146	USE caldyn_gcm_mod, ONLY : req_ps, req_mass
147	USE etat0_mod
148	USE guided_mod
149	USE advect_tracer_mod
150	USE physics_mod
151	USE mpipara
152	USE omp_para
153	USE trace
154	USE transfert_mod
155	USE check_conserve_mod
156	USE xios_mod
157	USE output_field_mod
158	USE write_etat0_mod
159	IMPLICIT NONE
160	REAL(rstd),POINTER :: q(:,:,:)
161	REAL(rstd),POINTER :: phis(:), ps(:) ,psm1(:), psm2(:), dps(:)
162	REAL(rstd),POINTER :: u(:,:) , um1(:,:), um2(:,:), du(:,:)
163	REAL(rstd),POINTER :: rhodz(:,:), mass(:,:), massm1(:,:), massm2(:,:), dmass(:,:)
164	REAL(rstd),POINTER :: theta_rhodz(:,:) , theta_rhodzm1(:,:), theta_rhodzm2(:,:), dtheta_rhodz(:,:)
165	REAL(rstd),POINTER :: hflux(:,:),wflux(:,:),hfluxt(:,:),wfluxt(:,:)
166
167	INTEGER :: ind
168	INTEGER :: it,i,j,n, stage
169	LOGICAL :: fluxt_zero(ndomain) ! set to .TRUE. to start accumulating fluxes in time
170	LOGICAL, PARAMETER :: check=.FALSE.
171	INTEGER :: start_clock
172	INTEGER :: stop_clock
173	INTEGER :: rate_clock
174
175
176	! CALL write_etat0(f_ps, f_phis,f_theta_rhodz,f_u,f_q)
177	! CALL read_start(f_ps,f_mass,f_phis,f_theta_rhodz,f_u,f_q)
178	! CALL write_restart(f_ps,f_mass,f_phis,f_theta_rhodz,f_u,f_q)
179
180	CALL caldyn_BC(f_phis, f_wflux) ! set constant values in first/last interfaces
181
182	!!$OMP BARRIER
183	DO ind=1,ndomain
184	IF (.NOT. assigned_domain(ind)) CYCLE
185	CALL swap_dimensions(ind)
186	CALL swap_geometry(ind)
187	rhodz=f_rhodz(ind); mass=f_mass(ind); ps=f_ps(ind)
188	IF(caldyn_eta==eta_mass) THEN
189	CALL compute_rhodz(.TRUE., ps, rhodz) ! save rhodz for transport scheme before dynamics update ps
190	ELSE
191	rhodz(:,:)=mass(:,:)
192	END IF
193	END DO
194	fluxt_zero=.TRUE.
195
196	!$OMP MASTER
197	CALL SYSTEM_CLOCK(start_clock)
198	!$OMP END MASTER
199
200	CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,itau0)
201
202	CALL trace_on
203
204	DO it=itau0+1,itau0+itaumax
205
206	IF (xios_output) CALL xios_update_calendar(it)
207	IF (it==itau0+1 .OR. MOD(it,itau_sync)==0) THEN
208	CALL send_message(f_ps,req_ps0)
209	CALL wait_message(req_ps0)
210	CALL send_message(f_mass,req_mass0)
211	CALL wait_message(req_mass0)
212	CALL send_message(f_theta_rhodz,req_theta_rhodz0)
213	CALL wait_message(req_theta_rhodz0)
214	CALL send_message(f_u,req_u0)
215	CALL wait_message(req_u0)
216	CALL send_message(f_q,req_q0)
217	CALL wait_message(req_q0)
218
219	! CALL wait_message(req_ps0)
220	! CALL wait_message(req_mass0)
221	! CALL wait_message(req_theta_rhodz0)
222	! CALL wait_message(req_u0)
223	! CALL wait_message(req_q0)
224	ENDIF
225
226	!$OMP MASTER
227	IF (is_mpi_root) PRINT ,"It No :",It," t :",dtIt
228	!$OMP END MASTER
229	IF (mod(it,itau_out)==0 ) THEN
230	CALL update_time_counter(dt*it)
231	CALL output_field("q",f_q)
232	CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it)
233	ENDIF
234
235	CALL guided(it*dt,f_ps,f_theta_rhodz,f_u,f_q)
236
237	DO stage=1,nb_stage
238	CALL caldyn((stage==1) .AND. (MOD(it,itau_out)==0), &
239	f_phis,f_ps,f_mass,f_theta_rhodz,f_u, f_q, &
240	f_hflux, f_wflux, f_dps, f_dmass, f_dtheta_rhodz, f_du)
241	SELECT CASE (scheme)
242	CASE(euler)
243	CALL euler_scheme(.TRUE.)
244	CASE (rk4)
245	CALL rk_scheme(stage)
246	CASE (mlf)
247	CALL leapfrog_matsuno_scheme(stage)
248	CASE DEFAULT
249	STOP
250	END SELECT
251	END DO
252
253	IF (MOD(it,itau_dissip)==0) THEN
254	! CALL send_message(f_ps,req_ps)
255	! CALL wait_message(req_ps)
256
257	IF(caldyn_eta==eta_mass) THEN
258	DO ind=1,ndomain
259	IF (.NOT. assigned_domain(ind)) CYCLE
260	CALL swap_dimensions(ind)
261	CALL swap_geometry(ind)
262	mass=f_mass(ind); ps=f_ps(ind);
263	CALL compute_rhodz(.TRUE., ps, mass)
264	END DO
265	ENDIF
266	! CALL send_message(f_mass,req_mass)
267	! CALL wait_message(req_mass)
268	CALL dissip(f_u,f_du,f_mass,f_phis, f_theta_rhodz,f_dtheta_rhodz)
269	! CALL send_message(f_mass,req_mass)
270	! CALL wait_message(req_mass)
271	CALL euler_scheme(.FALSE.) ! update only u, theta
272	END IF
273
274	IF(MOD(it,itau_adv)==0) THEN
275
276	CALL advect_tracer(f_hfluxt,f_wfluxt,f_u, f_q,f_rhodz) ! update q and rhodz after RK step
277	fluxt_zero=.TRUE.
278
279	! FIXME : check that rhodz is consistent with ps
280	IF (check) THEN
281	DO ind=1,ndomain
282	IF (.NOT. assigned_domain(ind)) CYCLE
283	CALL swap_dimensions(ind)
284	CALL swap_geometry(ind)
285	rhodz=f_rhodz(ind); ps=f_ps(ind);
286	CALL compute_rhodz(.FALSE., ps, rhodz)
287	END DO
288	ENDIF
289
290	END IF
291
292
293	IF (MOD(it,itau_physics)==0) THEN
294	CALL physics(it,f_phis, f_ps, f_theta_rhodz, f_u, f_q)
295	ENDIF
296
297	ENDDO
298
299	CALL write_etat0(itau0+itaumax,f_ps, f_phis,f_theta_rhodz,f_u,f_q)
300
301	CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it)
302
303	!$OMP MASTER
304	CALL SYSTEM_CLOCK(stop_clock)
305	CALL SYSTEM_CLOCK(count_rate=rate_clock)
306
307	IF (mpi_rank==0) THEN
308	PRINT ,"Time elapsed : ",(stop_clock-start_clock)1./rate_clock
309	ENDIF
310	!$OMP END MASTER
311
312	CONTAINS
313
314	SUBROUTINE Euler_scheme(with_dps)
315	IMPLICIT NONE
316	LOGICAL :: with_dps
317	INTEGER :: ind
318	INTEGER :: i,j,ij,l
319	CALL trace_start("Euler_scheme")
320
321	DO ind=1,ndomain
322	IF (.NOT. assigned_domain(ind)) CYCLE
323	CALL swap_dimensions(ind)
324	CALL swap_geometry(ind)
325
326	IF(with_dps) THEN ! update ps/mass
327	IF(caldyn_eta==eta_mass) THEN ! update ps
328	ps=f_ps(ind) ; dps=f_dps(ind) ;
329	IF (omp_first) THEN
330	!$SIMD
331	DO ij=ij_begin,ij_end
332	ps(ij)=ps(ij)+dt*dps(ij)
333	ENDDO
334	ENDIF
335	ELSE ! update mass
336	mass=f_mass(ind) ; dmass=f_dmass(ind) ;
337	DO l=1,llm
338	!$SIMD
339	DO ij=ij_begin,ij_end
340	mass(ij,l)=mass(ij,l)+dt*dmass(ij,l)
341	ENDDO
342	END DO
343	END IF
344
345	hflux=f_hflux(ind); hfluxt=f_hfluxt(ind)
346	wflux=f_wflux(ind); wfluxt=f_wfluxt(ind)
347	CALL accumulate_fluxes(hflux,wflux,hfluxt,wfluxt,dt,fluxt_zero(ind))
348	END IF ! update ps/mass
349
350	u=f_u(ind) ; theta_rhodz=f_theta_rhodz(ind)
351	du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind)
352
353	DO l=ll_begin,ll_end
354	!$SIMD
355	DO ij=ij_begin,ij_end
356	u(ij+u_right,l)=u(ij+u_right,l)+dt*du(ij+u_right,l)
357	u(ij+u_lup,l)=u(ij+u_lup,l)+dt*du(ij+u_lup,l)
358	u(ij+u_ldown,l)=u(ij+u_ldown,l)+dt*du(ij+u_ldown,l)
359	theta_rhodz(ij,l)=theta_rhodz(ij,l)+dt*dtheta_rhodz(ij,l)
360	ENDDO
361	ENDDO
362	ENDDO
363
364	CALL trace_end("Euler_scheme")
365
366	END SUBROUTINE Euler_scheme
367
368	SUBROUTINE RK_scheme(stage)
369	IMPLICIT NONE
370	INTEGER :: ind, stage
371	REAL(rstd), DIMENSION(4), PARAMETER :: coef = (/ .25, 1./3., .5, 1. /)
372	REAL(rstd) :: tau
373	INTEGER :: i,j,ij,l
374
375	CALL trace_start("RK_scheme")
376
377	tau = dt*coef(stage)
378
379	! if mass coordinate, deal with ps first on one core
380	IF(caldyn_eta==eta_mass) THEN
381	IF (omp_first) THEN
382
383	DO ind=1,ndomain
384	IF (.NOT. assigned_domain(ind)) CYCLE
385	CALL swap_dimensions(ind)
386	CALL swap_geometry(ind)
387	ps=f_ps(ind) ; psm1=f_psm1(ind) ; dps=f_dps(ind)
388
389	IF (stage==1) THEN ! first stage : save model state in XXm1
390	!$SIMD
391	DO ij=ij_begin,ij_end
392	psm1(ij)=ps(ij)
393	ENDDO
394	ENDIF
395
396	! updates are of the form x1 := x0 + tau*f(x1)
397	!$SIMD
398	DO ij=ij_begin,ij_end
399	ps(ij)=psm1(ij)+tau*dps(ij)
400	ENDDO
401	ENDDO
402	ENDIF
403	! CALL send_message(f_ps,req_ps)
404	!ym no overlap for now
405	! CALL wait_message(req_ps)
406
407	ELSE ! Lagrangian coordinate, deal with mass
408	DO ind=1,ndomain
409	IF (.NOT. assigned_domain(ind)) CYCLE
410	CALL swap_dimensions(ind)
411	CALL swap_geometry(ind)
412	mass=f_mass(ind); dmass=f_dmass(ind); massm1=f_massm1(ind)
413
414	IF (stage==1) THEN ! first stage : save model state in XXm1
415	DO l=ll_begin,ll_end
416	!$SIMD
417	DO ij=ij_begin,ij_end
418	massm1(ij,l)=mass(ij,l)
419	ENDDO
420	ENDDO
421	END IF
422
423	! updates are of the form x1 := x0 + tau*f(x1)
424	DO l=ll_begin,ll_end
425	!$SIMD
426	DO ij=ij_begin,ij_end
427	mass(ij,l)=massm1(ij,l)+tau*dmass(ij,l)
428	ENDDO
429	ENDDO
430	END DO
431	! CALL send_message(f_mass,req_mass)
432	!ym no overlap for now
433	! CALL wait_message(req_mass)
434
435	END IF
436
437	! now deal with other prognostic variables
438	DO ind=1,ndomain
439	IF (.NOT. assigned_domain(ind)) CYCLE
440	CALL swap_dimensions(ind)
441	CALL swap_geometry(ind)
442	u=f_u(ind) ; du=f_du(ind) ; um1=f_um1(ind)
443	theta_rhodz=f_theta_rhodz(ind)
444	theta_rhodzm1=f_theta_rhodzm1(ind)
445	dtheta_rhodz=f_dtheta_rhodz(ind)
446
447	IF (stage==1) THEN ! first stage : save model state in XXm1
448	DO l=ll_begin,ll_end
449	!$SIMD
450	DO ij=ij_begin,ij_end
451	um1(ij+u_right,l)=u(ij+u_right,l)
452	um1(ij+u_lup,l)=u(ij+u_lup,l)
453	um1(ij+u_ldown,l)=u(ij+u_ldown,l)
454	theta_rhodzm1(ij,l)=theta_rhodz(ij,l)
455	ENDDO
456	ENDDO
457	END IF
458
459	DO l=ll_begin,ll_end
460	!$SIMD
461	DO ij=ij_begin,ij_end
462	u(ij+u_right,l)=um1(ij+u_right,l)+tau*du(ij+u_right,l)
463	u(ij+u_lup,l)=um1(ij+u_lup,l)+tau*du(ij+u_lup,l)
464	u(ij+u_ldown,l)=um1(ij+u_ldown,l)+tau*du(ij+u_ldown,l)
465	theta_rhodz(ij,l)=theta_rhodzm1(ij,l)+tau*dtheta_rhodz(ij,l)
466	ENDDO
467	ENDDO
468
469	IF(stage==nb_stage) THEN ! accumulate mass fluxes at last stage
470	hflux=f_hflux(ind); hfluxt=f_hfluxt(ind)
471	wflux=f_wflux(ind); wfluxt=f_wfluxt(ind)
472	CALL accumulate_fluxes(hflux,wflux, hfluxt,wfluxt, dt,fluxt_zero(ind))
473	END IF
474	END DO
475
476	CALL trace_end("RK_scheme")
477
478	END SUBROUTINE RK_scheme
479
480	SUBROUTINE leapfrog_matsuno_scheme(stage)
481	IMPLICIT NONE
482	INTEGER :: ind, stage
483	REAL :: tau
484
485	CALL trace_start("leapfrog_matsuno_scheme")
486
487	tau = dt/nb_stage
488	DO ind=1,ndomain
489	IF (.NOT. assigned_domain(ind)) CYCLE
490	CALL swap_dimensions(ind)
491	CALL swap_geometry(ind)
492
493	ps=f_ps(ind) ; u=f_u(ind) ; theta_rhodz=f_theta_rhodz(ind)
494	psm1=f_psm1(ind) ; um1=f_um1(ind) ; theta_rhodzm1=f_theta_rhodzm1(ind)
495	psm2=f_psm2(ind) ; um2=f_um2(ind) ; theta_rhodzm2=f_theta_rhodzm2(ind)
496	dps=f_dps(ind) ; du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind)
497
498
499	IF (stage==1) THEN
500	psm1(:)=ps(:) ; um1(:,:)=u(:,:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:)
501	ps(:)=psm1(:)+tau*dps(:)
502	u(:,:)=um1(:,:)+tau*du(:,:)
503	theta_rhodz(:,:)=theta_rhodzm1(:,:)+tau*dtheta_rhodz(:,:)
504
505	ELSE IF (stage==2) THEN
506
507	ps(:)=psm1(:)+tau*dps(:)
508	u(:,:)=um1(:,:)+tau*du(:,:)
509	theta_rhodz(:,:)=theta_rhodzm1(:,:)+tau*dtheta_rhodz(:,:)
510
511	psm2(:)=psm1(:) ; theta_rhodzm2(:,:)=theta_rhodzm1(:,:) ; um2(:,:)=um1(:,:)
512	psm1(:)=ps(:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:) ; um1(:,:)=u(:,:)
513
514	ELSE
515
516	ps(:)=psm2(:)+2taudps(:)
517	u(:,:)=um2(:,:)+2taudu(:,:)
518	theta_rhodz(:,:)=theta_rhodzm2(:,:)+2taudtheta_rhodz(:,:)
519
520	psm2(:)=psm1(:) ; theta_rhodzm2(:,:)=theta_rhodzm1(:,:) ; um2(:,:)=um1(:,:)
521	psm1(:)=ps(:) ; theta_rhodzm1(:,:)=theta_rhodz(:,:) ; um1(:,:)=u(:,:)
522
523	ENDIF
524
525	ENDDO
526	CALL trace_end("leapfrog_matsuno_scheme")
527
528	END SUBROUTINE leapfrog_matsuno_scheme
529
530	END SUBROUTINE timeloop
531
532	SUBROUTINE accumulate_fluxes(hflux,wflux, hfluxt,wfluxt, tau,fluxt_zero)
533	USE icosa
534	USE omp_para
535	USE disvert_mod
536	IMPLICIT NONE
537	REAL(rstd), INTENT(IN) :: hflux(3iimjjm,llm), wflux(iim*jjm,llm+1)
538	REAL(rstd), INTENT(INOUT) :: hfluxt(3iimjjm,llm), wfluxt(iim*jjm,llm+1)
539	REAL(rstd), INTENT(IN) :: tau
540	LOGICAL, INTENT(INOUT) :: fluxt_zero
541	INTEGER :: l,i,j,ij
542
543	IF(fluxt_zero) THEN
544
545	fluxt_zero=.FALSE.
546
547	DO l=ll_begin,ll_end
548	!$SIMD
549	DO ij=ij_begin_ext,ij_end_ext
550	hfluxt(ij+u_right,l) = tau*hflux(ij+u_right,l)
551	hfluxt(ij+u_lup,l) = tau*hflux(ij+u_lup,l)
552	hfluxt(ij+u_ldown,l) = tau*hflux(ij+u_ldown,l)
553	ENDDO
554	ENDDO
555
556	IF(caldyn_eta==eta_mass) THEN ! no need for vertical fluxes if eta_lag
557	DO l=ll_begin,ll_endp1
558	!$SIMD
559	DO ij=ij_begin,ij_end
560	wfluxt(ij,l) = tau*wflux(ij,l)
561	ENDDO
562	ENDDO
563	END IF
564
565	ELSE
566
567	DO l=ll_begin,ll_end
568	!$SIMD
569	DO ij=ij_begin_ext,ij_end_ext
570	hfluxt(ij+u_right,l) = hfluxt(ij+u_right,l)+tau*hflux(ij+u_right,l)
571	hfluxt(ij+u_lup,l) = hfluxt(ij+u_lup,l)+tau*hflux(ij+u_lup,l)
572	hfluxt(ij+u_ldown,l) = hfluxt(ij+u_ldown,l)+tau*hflux(ij+u_ldown,l)
573	ENDDO
574	ENDDO
575
576	IF(caldyn_eta==eta_mass) THEN ! no need for vertical fluxes if eta_lag
577	DO l=ll_begin,ll_endp1
578	!$SIMD
579	DO ij=ij_begin,ij_end
580	wfluxt(ij,l) = wfluxt(ij,l)+tau*wflux(ij,l)
581	ENDDO
582	ENDDO
583	END IF
584
585	END IF
586
587	END SUBROUTINE accumulate_fluxes
588
589	! FUNCTION maxval_i(p)
590	! USE icosa
591	! IMPLICIT NONE
592	! REAL(rstd), DIMENSION(iim*jjm) :: p
593	! REAL(rstd) :: maxval_i
594	! INTEGER :: j, ij
595	!
596	! maxval_i=p((jj_begin-1)*iim+ii_begin)
597	!
598	! DO j=jj_begin-1,jj_end+1
599	! ij=(j-1)*iim
600	! maxval_i = MAX(maxval_i, MAXVAL(p(ij+ii_begin:ij+ii_end)))
601	! END DO
602	! END FUNCTION maxval_i
603
604	! FUNCTION maxval_ik(p)
605	! USE icosa
606	! IMPLICIT NONE
607	! REAL(rstd) :: p(iim*jjm, llm)
608	! REAL(rstd) :: maxval_ik(llm)
609	! INTEGER :: l,j, ij
610	!
611	! DO l=1,llm
612	! maxval_ik(l)=p((jj_begin-1)*iim+ii_begin,l)
613	! DO j=jj_begin-1,jj_end+1
614	! ij=(j-1)*iim
615	! maxval_ik(l) = MAX(maxval_ik(l), MAXVAL(p(ij+ii_begin:ij+ii_end,l)))
616	! END DO
617	! END DO
618	! END FUNCTION maxval_ik
619
620	END MODULE timeloop_gcm_mod

Note: See TracBrowser for help on using the repository browser.

Download in other formats: