Context Navigation

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

Last change on this file since 178 was 174, checked in by ymipsl, 11 years ago

Transform 2 loops on i and j in one loop ij for efficiency, vectorization and future GPU programing

YM

File size: 18.9 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 :: itau_sync=10
10
11	TYPE(t_message) :: req_ps0, req_mass0, req_theta_rhodz0, req_u0, req_q0
12
13	TYPE(t_field),POINTER :: f_q(:)
14	TYPE(t_field),POINTER :: f_rhodz(:), f_mass(:), f_massm1(:), f_massm2(:), f_dmass(:)
15	TYPE(t_field),POINTER :: f_phis(:), f_ps(:),f_psm1(:), f_psm2(:), f_dps(:)
16	TYPE(t_field),POINTER :: f_u(:),f_um1(:),f_um2(:), f_du(:)
17	TYPE(t_field),POINTER :: f_theta_rhodz(:),f_theta_rhodzm1(:),f_theta_rhodzm2(:), f_dtheta_rhodz(:)
18	TYPE(t_field),POINTER :: f_hflux(:), f_wflux(:), f_hfluxt(:), f_wfluxt(:)
19
20	INTEGER :: nb_stage, matsuno_period, scheme
21
22	REAL(rstd),SAVE :: jD_cur, jH_cur
23	REAL(rstd),SAVE :: start_time
24
25	CONTAINS
26
27	SUBROUTINE init_timeloop
28	USE icosa
29	USE dissip_gcm_mod
30	USE caldyn_mod
31	USE etat0_mod
32	USE disvert_mod
33	USE guided_mod
34	USE advect_tracer_mod
35	USE physics_mod
36	USE mpipara
37	USE omp_para
38	USE trace
39	USE transfert_mod
40	USE check_conserve_mod
41	USE ioipsl
42	USE output_field_mod
43	USE write_field
44	IMPLICIT NONE
45
46	CHARACTER(len=255) :: def
47
48	!----------------------------------------------------
49	IF (TRIM(time_style)=='lmd') Then
50
51	day_step=180
52	CALL getin('day_step',day_step)
53
54	ndays=1
55	CALL getin('ndays',ndays)
56
57	dt = daysec/REAL(day_step)
58	itaumax = ndays*day_step
59
60	calend = 'earth_360d'
61	CALL getin('calend', calend)
62
63	day_ini = 0
64	CALL getin('day_ini',day_ini)
65
66	day_end = 0
67	CALL getin('day_end',day_end)
68
69	annee_ref = 1998
70	CALL getin('annee_ref',annee_ref)
71
72	start_time = 0
73	CALL getin('start_time',start_time)
74
75
76	write_period=0
77	CALL getin('write_period',write_period)
78
79	write_period=write_period/scale_factor
80	itau_out=FLOOR(write_period/dt)
81
82	PRINT *, 'Output frequency (scaled) set to ',write_period, ' : itau_out = ',itau_out
83
84	mois = 1 ; heure = 0.
85	call ymds2ju(annee_ref, mois, day_ref, heure, jD_ref)
86	jH_ref = jD_ref - int(jD_ref)
87	jD_ref = int(jD_ref)
88
89	CALL ioconf_startdate(INT(jD_ref),jH_ref)
90	write(,)'annee_ref, mois, day_ref, heure, jD_ref'
91	write(,)annee_ref, mois, day_ref, heure, jD_ref
92	write(,)"ndays,day_step,itaumax,dt======>"
93	write(,)ndays,day_step,itaumax,dt
94	call ju2ymds(jD_ref+jH_ref,an, mois, jour, heure)
95	write(,)'jD_ref+jH_ref,an, mois, jour, heure'
96	write(,)jD_ref+jH_ref,an, mois, jour, heure
97	day_end = day_ini + ndays
98	END IF
99	!----------------------------------------------------
100
101	IF (xios_output) itau_out=1
102
103	! Time-independant orography
104	CALL allocate_field(f_phis,field_t,type_real,name='phis')
105	! Trends
106	CALL allocate_field(f_du,field_u,type_real,llm,name='du')
107	CALL allocate_field(f_dtheta_rhodz,field_t,type_real,llm,name='dtheta_rhodz')
108	! Model state at current time step (RK/MLF/Euler)
109	CALL allocate_field(f_ps,field_t,type_real, name='ps')
110	CALL allocate_field(f_mass,field_t,type_real,llm,name='mass')
111	CALL allocate_field(f_u,field_u,type_real,llm,name='u')
112	CALL allocate_field(f_theta_rhodz,field_t,type_real,llm,name='theta_rhodz')
113	! Model state at previous time step (RK/MLF)
114	CALL allocate_field(f_um1,field_u,type_real,llm,name='um1')
115	CALL allocate_field(f_theta_rhodzm1,field_t,type_real,llm,name='theta_rhodzm1')
116	! Tracers
117	CALL allocate_field(f_q,field_t,type_real,llm,nqtot)
118	CALL allocate_field(f_rhodz,field_t,type_real,llm,name='rhodz')
119	! Mass fluxes
120	CALL allocate_field(f_hflux,field_u,type_real,llm) ! instantaneous mass fluxes
121	CALL allocate_field(f_hfluxt,field_u,type_real,llm) ! mass "fluxes" accumulated in time
122	CALL allocate_field(f_wflux,field_t,type_real,llm+1) ! vertical mass fluxes
123	CALL allocate_field(f_dmass,field_t,type_real,llm, name='dmass')
124
125	IF(caldyn_eta == eta_mass) THEN ! eta = mass coordinate (default)
126	CALL allocate_field(f_dps,field_t,type_real,name='dps')
127	CALL allocate_field(f_psm1,field_t,type_real,name='psm1')
128	CALL allocate_field(f_wfluxt,field_t,type_real,llm+1,name='wfluxt')
129	! the following are unused but must point to something
130	! f_massm1 => f_mass
131	ELSE ! eta = Lagrangian vertical coordinate
132	CALL allocate_field(f_massm1,field_t,type_real,llm, name='massm1')
133	! the following are unused but must point to something
134	f_wfluxt => f_wflux
135	f_dps => f_phis
136	f_psm1 => f_phis
137	END IF
138
139	def='runge_kutta'
140	CALL getin('scheme',def)
141
142	SELECT CASE (TRIM(def))
143	CASE('euler')
144	scheme=euler
145	nb_stage=1
146	CASE ('runge_kutta')
147	scheme=rk4
148	nb_stage=4
149	CASE ('leapfrog_matsuno')
150	scheme=mlf
151	matsuno_period=5
152	CALL getin('matsuno_period',matsuno_period)
153	nb_stage=matsuno_period+1
154	! Model state 2 time steps ago (MLF)
155	CALL allocate_field(f_theta_rhodzm2,field_t,type_real,llm)
156	CALL allocate_field(f_um2,field_u,type_real,llm)
157	IF(caldyn_eta == eta_mass) THEN ! eta = mass coordinate (default)
158	CALL allocate_field(f_psm2,field_t,type_real)
159	! the following are unused but must point to something
160	f_massm2 => f_mass
161	ELSE ! eta = Lagrangian vertical coordinate
162	CALL allocate_field(f_massm2,field_t,type_real,llm)
163	! the following are unused but must point to something
164	f_psm2 => f_phis
165	END IF
166
167	CASE default
168	PRINT*,'Bad selector for variable scheme : <', TRIM(def), &
169	' > options are <euler>, <runge_kutta>, <leapfrog_matsuno>'
170	STOP
171	END SELECT
172
173
174	CALL init_dissip
175	CALL init_caldyn
176	CALL init_guided
177	CALL init_advect_tracer
178	CALL init_check_conserve
179	CALL init_physics
180	CALL etat0(f_ps,f_mass,f_phis,f_theta_rhodz,f_u, f_q)
181
182	CALL transfert_request(f_phis,req_i0)
183	CALL transfert_request(f_phis,req_i1)
184	CALL writefield("phis",f_phis,once=.TRUE.)
185
186	CALL init_message(f_ps,req_i0,req_ps0)
187	CALL init_message(f_mass,req_i0,req_mass0)
188	CALL init_message(f_theta_rhodz,req_i0,req_theta_rhodz0)
189	CALL init_message(f_u,req_e0_vect,req_u0)
190	CALL init_message(f_q,req_i0,req_q0)
191
192	END SUBROUTINE init_timeloop
193
194	SUBROUTINE timeloop
195	USE icosa
196	USE dissip_gcm_mod
197	USE disvert_mod
198	USE caldyn_mod
199	USE caldyn_gcm_mod, ONLY : req_ps, req_mass
200	USE etat0_mod
201	USE guided_mod
202	USE advect_tracer_mod
203	USE physics_mod
204	USE mpipara
205	USE omp_para
206	USE trace
207	USE transfert_mod
208	USE check_conserve_mod
209	USE xios_mod
210	USE output_field_mod
211	IMPLICIT NONE
212	REAL(rstd),POINTER :: q(:,:,:)
213	REAL(rstd),POINTER :: phis(:), ps(:) ,psm1(:), psm2(:), dps(:)
214	REAL(rstd),POINTER :: u(:,:) , um1(:,:), um2(:,:), du(:,:)
215	REAL(rstd),POINTER :: rhodz(:,:), mass(:,:), massm1(:,:), massm2(:,:), dmass(:,:)
216	REAL(rstd),POINTER :: theta_rhodz(:,:) , theta_rhodzm1(:,:), theta_rhodzm2(:,:), dtheta_rhodz(:,:)
217	REAL(rstd),POINTER :: hflux(:,:),wflux(:,:),hfluxt(:,:),wfluxt(:,:)
218
219	INTEGER :: ind
220	INTEGER :: it,i,j,n, stage
221	LOGICAL :: fluxt_zero(ndomain) ! set to .TRUE. to start accumulating fluxes in time
222	LOGICAL, PARAMETER :: check=.FALSE.
223
224	CALL caldyn_BC(f_phis, f_wflux) ! set constant values in first/last interfaces
225
226	!$OMP BARRIER
227	DO ind=1,ndomain
228	CALL swap_dimensions(ind)
229	CALL swap_geometry(ind)
230	rhodz=f_rhodz(ind); mass=f_mass(ind); ps=f_ps(ind)
231	IF(caldyn_eta==eta_mass) THEN
232	CALL compute_rhodz(.TRUE., ps, rhodz) ! save rhodz for transport scheme before dynamics update ps
233	ELSE
234	rhodz(:,:)=mass(:,:)
235	END IF
236	END DO
237	fluxt_zero=.TRUE.
238
239	DO it=0,itaumax
240
241	CALL xios_update_calendar(it)
242	IF (MOD(it,itau_sync)==0) THEN
243	CALL send_message(f_ps,req_ps0)
244	CALL send_message(f_mass,req_mass0)
245	CALL send_message(f_theta_rhodz,req_theta_rhodz0)
246	CALL send_message(f_u,req_u0)
247	CALL send_message(f_q,req_q0)
248	CALL wait_message(req_ps0)
249	CALL wait_message(req_mass0)
250	CALL wait_message(req_theta_rhodz0)
251	CALL wait_message(req_u0)
252	CALL wait_message(req_q0)
253	ENDIF
254
255	! IF (is_mpi_root) PRINT ,"It No :",It," t :",dtIt
256	IF (mod(it,itau_out)==0 ) THEN
257	CALL update_time_counter(dt*it)
258	CALL output_field("q",f_q)
259	CALL check_conserve(f_ps,f_dps,f_u,f_theta_rhodz,f_phis,it)
260	ENDIF
261
262	CALL guided(it*dt,f_ps,f_theta_rhodz,f_u,f_q)
263
264	DO stage=1,nb_stage
265	CALL caldyn((stage==1) .AND. (MOD(it,itau_out)==0), &
266	f_phis,f_ps,f_mass,f_theta_rhodz,f_u, f_q, &
267	f_hflux, f_wflux, f_dps, f_dmass, f_dtheta_rhodz, f_du)
268	SELECT CASE (scheme)
269	CASE(euler)
270	CALL euler_scheme(.TRUE.)
271	CASE (rk4)
272	CALL rk_scheme(stage)
273	CASE (mlf)
274	CALL leapfrog_matsuno_scheme(stage)
275	CASE DEFAULT
276	STOP
277	END SELECT
278	END DO
279
280	IF (MOD(it+1,itau_dissip)==0) THEN
281	IF(caldyn_eta==eta_mass) THEN
282	DO ind=1,ndomain
283	CALL swap_dimensions(ind)
284	CALL swap_geometry(ind)
285	mass=f_mass(ind); ps=f_ps(ind);
286	CALL compute_rhodz(.TRUE., ps, mass)
287	END DO
288	ENDIF
289	CALL dissip(f_u,f_du,f_mass,f_phis, f_theta_rhodz,f_dtheta_rhodz)
290	CALL euler_scheme(.FALSE.) ! update only u, theta
291	END IF
292
293	IF(MOD(it+1,itau_adv)==0) THEN
294
295	CALL advect_tracer(f_hfluxt,f_wfluxt,f_u, f_q,f_rhodz) ! update q and rhodz after RK step
296	fluxt_zero=.TRUE.
297
298	! FIXME : check that rhodz is consistent with ps
299	IF (check) THEN
300	DO ind=1,ndomain
301	CALL swap_dimensions(ind)
302	CALL swap_geometry(ind)
303	rhodz=f_rhodz(ind); ps=f_ps(ind);
304	CALL compute_rhodz(.FALSE., ps, rhodz)
305	END DO
306	ENDIF
307
308	END IF
309
310
311
312	!----------------------------------------------------
313	! jD_cur = jD_ref + day_ini - day_ref + it/day_step
314	! jH_cur = jH_ref + start_time + mod(it,day_step)/float(day_step)
315	! jD_cur = jD_cur + int(jH_cur)
316	! jH_cur = jH_cur - int(jH_cur)
317	CALL physics(it,jD_cur,jH_cur,f_phis, f_ps, f_theta_rhodz, f_u, f_q)
318	ENDDO
319
320
321	CONTAINS
322
323	SUBROUTINE Euler_scheme(with_dps)
324	IMPLICIT NONE
325	LOGICAL :: with_dps
326	INTEGER :: ind
327	INTEGER :: i,j,ij,l
328	CALL trace_start("Euler_scheme")
329
330	DO ind=1,ndomain
331	CALL swap_dimensions(ind)
332	CALL swap_geometry(ind)
333
334	IF(with_dps) THEN ! update ps/mass
335	IF(caldyn_eta==eta_mass) THEN ! update ps
336	ps=f_ps(ind) ; dps=f_dps(ind) ;
337	IF (omp_first) THEN
338	!$SIMD
339	DO ij=ij_begin,ij_end
340	ps(ij)=ps(ij)+dt*dps(ij)
341	ENDDO
342	ENDIF
343	ELSE ! update mass
344	mass=f_mass(ind) ; dmass=f_dmass(ind) ;
345	DO l=1,llm
346	!$SIMD
347	DO ij=ij_begin,ij_end
348	mass(ij,l)=mass(ij,l)+dt*dmass(ij,l)
349	ENDDO
350	END DO
351	END IF
352
353	hflux=f_hflux(ind); hfluxt=f_hfluxt(ind)
354	wflux=f_wflux(ind); wfluxt=f_wfluxt(ind)
355	CALL accumulate_fluxes(hflux,wflux,hfluxt,wfluxt,dt,fluxt_zero(ind))
356	END IF ! update ps/mass
357
358	u=f_u(ind) ; theta_rhodz=f_theta_rhodz(ind)
359	du=f_du(ind) ; dtheta_rhodz=f_dtheta_rhodz(ind)
360
361	DO l=ll_begin,ll_end
362	!$SIMD
363	DO ij=ij_begin,ij_end
364	u(ij+u_right,l)=u(ij+u_right,l)+dt*du(ij+u_right,l)
365	u(ij+u_lup,l)=u(ij+u_lup,l)+dt*du(ij+u_lup,l)
366	u(ij+u_ldown,l)=u(ij+u_ldown,l)+dt*du(ij+u_ldown,l)
367	theta_rhodz(ij,l)=theta_rhodz(ij,l)+dt*dtheta_rhodz(ij,l)
368	ENDDO
369	ENDDO
370	ENDDO
371
372	CALL trace_end("Euler_scheme")
373
374	END SUBROUTINE Euler_scheme
375
376	SUBROUTINE RK_scheme(stage)
377	IMPLICIT NONE
378	INTEGER :: ind, stage
379	REAL(rstd), DIMENSION(4), PARAMETER :: coef = (/ .25, 1./3., .5, 1. /)
380	REAL(rstd) :: tau
381	INTEGER :: i,j,ij,l
382
383	CALL trace_start("RK_scheme")
384
385	tau = dt*coef(stage)
386
387	! if mass coordinate, deal with ps first on one core
388	IF(caldyn_eta==eta_mass) THEN
389	IF (omp_first) THEN
390
391	DO ind=1,ndomain
392	CALL swap_dimensions(ind)
393	CALL swap_geometry(ind)
394	ps=f_ps(ind) ; psm1=f_psm1(ind) ; dps=f_dps(ind)
395
396	IF (stage==1) THEN ! first stage : save model state in XXm1
397	!$SIMD
398	DO ij=ij_begin,ij_end
399	psm1(ij)=ps(ij)
400	ENDDO
401	ENDIF
402
403	! updates are of the form x1 := x0 + tau*f(x1)
404	!$SIMD
405	DO ij=ij_begin,ij_end
406	ps(ij)=psm1(ij)+tau*dps(ij)
407	ENDDO
408	ENDDO
409	ENDIF
410	CALL send_message(f_ps,req_ps)
411
412	ELSE ! Lagrangian coordinate, deal with mass
413	DO ind=1,ndomain
414	CALL swap_dimensions(ind)
415	CALL swap_geometry(ind)
416	mass=f_mass(ind); dmass=f_dmass(ind); massm1=f_massm1(ind)
417
418	IF (stage==1) THEN ! first stage : save model state in XXm1
419	DO l=ll_begin,ll_end
420	!$SIMD
421	DO ij=ij_begin,ij_end
422	massm1(ij,l)=mass(ij,l)
423	ENDDO
424	ENDDO
425	END IF
426
427	! updates are of the form x1 := x0 + tau*f(x1)
428	DO l=ll_begin,ll_end
429	!$SIMD
430	DO ij=ij_begin,ij_end
431	mass(ij,l)=massm1(ij,l)+tau*dmass(ij,l)
432	ENDDO
433	ENDDO
434	END DO
435	CALL send_message(f_mass,req_mass)
436
437	END IF
438
439	! now deal with other prognostic variables
440	DO ind=1,ndomain
441	CALL swap_dimensions(ind)
442	CALL swap_geometry(ind)
443	u=f_u(ind) ; du=f_du(ind) ; um1=f_um1(ind)
444	theta_rhodz=f_theta_rhodz(ind)
445	theta_rhodzm1=f_theta_rhodzm1(ind)
446	dtheta_rhodz=f_dtheta_rhodz(ind)
447
448	IF (stage==1) THEN ! first stage : save model state in XXm1
449	DO l=ll_begin,ll_end
450	!$SIMD
451	DO ij=ij_begin,ij_end
452	um1(ij+u_right,l)=u(ij+u_right,l)
453	um1(ij+u_lup,l)=u(ij+u_lup,l)
454	um1(ij+u_ldown,l)=u(ij+u_ldown,l)
455	theta_rhodzm1(ij,l)=theta_rhodz(ij,l)
456	ENDDO
457	ENDDO
458	END IF
459
460	DO l=ll_begin,ll_end
461	!$SIMD
462	DO ij=ij_begin,ij_end
463	u(ij+u_right,l)=um1(ij+u_right,l)+tau*du(ij+u_right,l)
464	u(ij+u_lup,l)=um1(ij+u_lup,l)+tau*du(ij+u_lup,l)
465	u(ij+u_ldown,l)=um1(ij+u_ldown,l)+tau*du(ij+u_ldown,l)
466	theta_rhodz(ij,l)=theta_rhodzm1(ij,l)+tau*dtheta_rhodz(ij,l)
467	ENDDO
468	ENDDO
469
470	IF(stage==nb_stage) THEN ! accumulate mass fluxes at last stage
471	hflux=f_hflux(ind); hfluxt=f_hfluxt(ind)
472	wflux=f_wflux(ind); wfluxt=f_wfluxt(ind)
473	CALL accumulate_fluxes(hflux,wflux, hfluxt,wfluxt, dt,fluxt_zero(ind))
474	END IF
475	END DO
476
477	CALL trace_end("RK_scheme")
478
479	END SUBROUTINE RK_scheme
480
481	SUBROUTINE leapfrog_matsuno_scheme(stage)
482	IMPLICIT NONE
483	INTEGER :: ind, stage
484	REAL :: tau
485
486	CALL trace_start("leapfrog_matsuno_scheme")
487
488	tau = dt/nb_stage
489	DO ind=1,ndomain
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: