ROMS
Loading...
Searching...
No Matches
rp_main3d.F
Go to the documentation of this file.
1#include "cppdefs.h"
2#if defined TL_IOMS && defined SOLVE3D
3 SUBROUTINE rp_main3d (RunInterval)
4!
5!git $Id$
6!================================================== Hernan G. Arango ===
7! Copyright (c) 2002-2025 The ROMS Group Andrew M. Moore !
8! Licensed under a MIT/X style license !
9! See License_ROMS.md !
10!=======================================================================
11! !
12! This routine is the main driver for representers tangent linear !
13! ROMS when configure as a full 3D baroclinic ocean model. !
14! It advances forward the representer model equations for all !
15! nested grids, if any, for the specified time interval (seconds), !
16! RunInterval. !
17! !
18!=======================================================================
19!
20 USE mod_param
21 USE mod_parallel
22# if defined MODEL_COUPLING && defined MCT_LIB
23 USE mod_coupler
24# endif
25 USE mod_iounits
26 USE mod_scalars
27 USE mod_stepping
28!
29# ifdef ANA_VMIX
30 USE analytical_mod, ONLY : ana_vmix
31# endif
32 USE dateclock_mod, ONLY : time_string
33# ifdef TIDE_GENERATING_FORCES
34 USE equilibrium_tide_mod, ONLY : equilibrium_tide
35# endif
36# if defined ATM_COUPLING_NOT_YET && defined MCT_LIB
37 USE mct_coupler_mod, ONLY : ocn2atm_coupling
38# endif
39# if defined WAV_COUPLING_NOT_YET && defined MCT_LIB
40 USE mct_coupler_mod, ONLY : ocn2wav_coupling
41# endif
42# ifdef FORWARD_READ
43 USE omega_mod, ONLY : omega
44 USE set_depth_mod, ONLY : set_depth
46# endif
47 USE strings_mod, ONLY : founderror
48# ifdef BIOLOGY
49 USE rp_biology_mod, ONLY : rp_biology
50# endif
51# ifdef BBL_MODEL_NOT_YET
52!! USE rp_bbl_mod, ONLY : rp_bblm
53# endif
54# if defined BULK_FLUXES_NOT_YET && !defined PRIOR_BULK_FLUXES
56# endif
57# ifdef BVF_MIXING_NOT_YET
58!! USE rp_bvf_mix_mod, ONLY : rp_bvf_mix
59# endif
60 USE rp_diag_mod, ONLY : rp_diag
61# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
63# endif
64# ifdef GLS_MIXING_NOT_YET
65!! USE rp_gls_corstep_mod, ONLY : rp_gls_corstep
66!! USE rp_gls_prestep_mod, ONLY : rp_gls_prestep
67# endif
68# ifdef LMD_MIXING_NOT_YET
69!! USE rp_lmd_vmix_mod, ONLY : rp_lmd_vmix
70# endif
71# ifdef MY25_MIXING
72!! USE rp_my25_corstep_mod, ONLY : rp_my25_corstep
73!! USE rp_my25_prestep_mod, ONLY : rp_my25_prestep
74# endif
75# ifdef ADJUST_BOUNDARY
79# endif
80 USE rp_omega_mod, ONLY : rp_omega
82# ifdef NEARSHORE_MELLOR_NOT_YET
83!! USE rp_radiation_stress_mod, ONLY : rp_radiation_stress
84# endif
85# ifndef TS_FIXED
86 USE rp_rho_eos_mod, ONLY : rp_rho_eos
87# endif
88 USE rp_rhs3d_mod, ONLY : rp_rhs3d
89# ifdef SEDIMENT_NOT_YET
90!! USE rp_sediment_mod, ONLY : rp_sediment
91# endif
94# if defined SSH_TIDES_NOT_YET || defined UV_TIDES_NOT_YET
95!! USE rp_set_tides_mod, ONLY : rp_set_tides
96# endif
97 USE rp_set_vbc_mod, ONLY : rp_set_vbc
99 USE rp_step2d_mod, ONLY : rp_step2d
100# ifndef TS_FIXED
101 USE rp_step3d_t_mod, ONLY : rp_step3d_t
102# endif
104# ifdef FLOATS_NOT_YET
105!! USE rp_step_floats_mod, ONLY : rp_step_floats
106# endif
107# ifdef WEAK_CONSTRAINT
108 USE tl_forcing_mod, ONLY : tl_forcing
109# endif
110# ifdef RP_AVERAGES
111 USE tl_set_avg_mod, ONLY : tl_set_avg
112# endif
113!! USE wvelocity_mod, ONLY : wvelocity
114!
115 implicit none
116!
117! Imported variable declarations.
118!
119 real(dp), intent(in) :: RunInterval
120!
121! Local variable declarations.
122!
123 integer :: ng, tile
124 integer :: my_iif, next_indx1
125# ifdef FLOATS_NOT_YET
126 integer :: Lend, Lstr, chunk_size
127# endif
128!
129 real(r8) :: MaxDT, my_StepTime
130!
131 character (len=*), parameter :: MyFile = &
132 & __FILE__
133!
134!=======================================================================
135! Time-step representers tangent linear 3D primitive equations.
136!=======================================================================
137!
138 my_steptime=0.0_r8
139 maxdt=maxval(dt)
140
141 step_loop : DO WHILE (my_steptime.le.(runinterval+0.5_r8*maxdt))
142
143 my_steptime=my_steptime+maxdt
144!
145! Set time indices and time clock.
146!
147 DO ng=1,ngrids
148 iic(ng)=iic(ng)+1
149 nstp(ng)=1+mod(iic(ng)-ntstart(ng),2)
150 nnew(ng)=3-nstp(ng)
151 nrhs(ng)=nstp(ng)
152!$OMP MASTER
153 time(ng)=time(ng)+dt(ng)
154 tdays(ng)=time(ng)*sec2day
155 CALL time_string (time(ng), time_code(ng))
156!$OMP END MASTER
157 END DO
158!$OMP BARRIER
159!
160!-----------------------------------------------------------------------
161! Read in required data, if any, data from input NetCDF files.
162!-----------------------------------------------------------------------
163!
164 DO ng=1,ngrids
165!$OMP MASTER
166 CALL rp_get_data (ng)
167!$OMP END MASTER
168!$OMP BARRIER
169 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
170 END DO
171!
172!-----------------------------------------------------------------------
173! If applicable, process input data: time interpolate between data
174! snapshots. Compute BASIC STATE depths and thickness.
175!-----------------------------------------------------------------------
176!
177 DO ng=1,ngrids
178 DO tile=first_tile(ng),last_tile(ng),+1
179 CALL rp_set_data (ng, tile)
180# ifdef FORWARD_READ
181 CALL set_depth (ng, tile, irpm)
182# endif
183 END DO
184!$OMP BARRIER
185 END DO
186 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
187
188# ifdef FORWARD_READ
189!
190!-----------------------------------------------------------------------
191! Compute BASIC STATE horizontal mass fluxes (Hz*u/n and Hz*v/m).
192!-----------------------------------------------------------------------
193!
194 DO ng=1,ngrids
195 DO tile=last_tile(ng),first_tile(ng),-1
196 CALL set_massflux (ng, tile, irpm)
197 END DO
198!$OMP BARRIER
199 END DO
200# endif
201# ifdef WEAK_CONSTRAINT
202!
203!-----------------------------------------------------------------------
204! If appropriate, add convolved adjoint solution impulse forcing to
205! the representer model solution. Notice that the forcing is only
206! needed after finishing all inner loops. The forcing is continuous.
207! That is, it is time interpolated at every time-step from available
208! snapshots (FrequentImpulse=TRUE).
209!-----------------------------------------------------------------------
210!
211 DO ng=1,ngrids
212 IF (frequentimpulse(ng)) THEN
213 DO tile=first_tile(ng),last_tile(ng),+1
214 CALL tl_forcing (ng, tile, kstp(ng), nstp(ng))
215 END DO
216!$OMP BARRIER
217 END IF
218 END DO
219# endif
220!
221!-----------------------------------------------------------------------
222! On the first timestep, it computes the initial depths and level
223! thicknesses from the initial free-surface field. Additionally, it
224! initializes the tangent linear state variables for all time levels
225! and applies lateral boundary conditions.
226!-----------------------------------------------------------------------
227!
228 DO ng=1,ngrids
229 IF (iic(ng).eq.ntstart(ng)) THEN
230 CALL rp_post_initial (ng, irpm)
231 END IF
232 END DO
233!
234!-----------------------------------------------------------------------
235! Compute horizontal mass fluxes (Hz*u/n and Hz*v/m), density related
236! quatities and report global diagnostics. Compute BASIC STATE omega
237! vertical velocity.
238!-----------------------------------------------------------------------
239!
240 DO ng=1,ngrids
241 DO tile=first_tile(ng),last_tile(ng),+1
242 CALL rp_set_massflux (ng, tile, irpm)
243# ifndef TS_FIXED
244 CALL rp_rho_eos (ng, tile, irpm)
245# endif
246# ifdef TIDE_GENERATING_FORCES
247 CALL equilibrium_tide (ng, tile, itlm)
248# endif
249 CALL rp_diag (ng, tile)
250# ifdef FORWARD_READ
251 CALL omega (ng, tile, irpm)
252# endif
253 END DO
254!$OMP BARRIER
255 END DO
256 IF (founderror(exit_flag, noerror, __line__, myfile)) RETURN
257
258# if defined ATM_COUPLING_NOT_YET && defined MCT_LIB
259!
260!-----------------------------------------------------------------------
261! Couple to atmospheric model every CoupleSteps(Iatmos) timesteps: get
262! air/sea fluxes.
263!-----------------------------------------------------------------------
264!
265 DO ng=1,ngrids
266 IF ((iic(ng).ne.ntstart(ng)).and. &
267 & mod(iic(ng)-1,couplesteps(iatmos,ng)).eq.0) THEN
268 DO tile=last_tile(ng),first_tile(ng),-1
269 CALL ocn2atm_coupling (ng, tile)
270 END DO
271!$OMP BARRIER
272 END IF
273 END DO
274# endif
275
276# if defined WAV_COUPLING_NOT_YET && defined MCT_LIB
277!
278!-----------------------------------------------------------------------
279! Couple to waves model every CoupleSteps(Iwaves) timesteps: get
280! waves/sea fluxes.
281!-----------------------------------------------------------------------
282!
283 DO ng=1,ngrids
284 IF ((iic(ng).ne.ntstart(ng)).and. &
285 & mod(iic(ng),couplesteps(iwaves,ng)).eq.0) THEN
286 DO tile=first_tile(ng),last_tile(ng),+1
287 CALL ocn2wav_coupling (ng, tile)
288 END DO
289!$OMP BARRIER
290 END IF
291 END DO
292# endif
293
294# ifdef NEARSHORE_MELLOR_NOT_YET
295!
296!-----------------------------------------------------------------------
297! Compute radiation stress terms.
298!-----------------------------------------------------------------------
299!
300 DO ng=1,ngrids
301 DO tile=last_tile(ng),first_tile(ng),-1
302 CALL rp_radiation_stress (ng, tile)
303 END DO
304!$OMP BARRIER
305 END DO
306# endif
307!
308!-----------------------------------------------------------------------
309! Set fields for vertical boundary conditions. Process tidal forcing,
310! if any.
311!-----------------------------------------------------------------------
312!
313 DO ng=1,ngrids
314 DO tile=first_tile(ng),last_tile(ng),+1
315# if defined BULK_FLUXES_NOT_YET && !defined PRIOR_BULK_FLUXES
316 CALL rp_bulk_flux (ng, tile)
317# endif
318# ifdef BBL_MODEL_NOT_YET
319 CALL rp_bblm (ng, tile)
320# endif
321 CALL rp_set_vbc (ng, tile)
322# if defined SSH_TIDES_NOT_YET || defined UV_TIDES_NOT_YET
323 CALL rp_set_tides (ng, tile)
324# endif
325 END DO
326!$OMP BARRIER
327 END DO
328
329# ifdef ADJUST_BOUNDARY
330!
331!-----------------------------------------------------------------------
332! Interpolate open boundary increments and adjust open boundaries.
333! Skip the last output timestep.
334!-----------------------------------------------------------------------
335!
336 DO ng=1,ngrids
337 IF ((nrun.ne.1).and.(iic(ng).lt.(ntend(ng)+1))) THEN
338 DO tile=first_tile(ng),last_tile(ng),+1
339 CALL rp_obc_adjust (ng, tile, lbinp(ng))
340 CALL rp_set_depth_bry (ng, tile, irpm)
341 CALL rp_obc2d_adjust (ng, tile, lbinp(ng))
342 END DO
343!$OMP BARRIER
344 END IF
345 END DO
346# endif
347
348# if defined ADJUST_STFLUX || defined ADJUST_WSTRESS
349!
350!-----------------------------------------------------------------------
351! Interpolate surface forcing increments and adjust surface forcing.
352! Skip the last output timestep.
353!-----------------------------------------------------------------------
354!
355 DO ng=1,ngrids
356 IF (iic(ng).lt.(ntend(ng)+1)) THEN
357 DO tile=first_tile(ng),last_tile(ng),+1
358 CALL rp_frc_adjust (ng, tile, lfinp(ng))
359 END DO
360!$OMP BARRIER
361 END IF
362 END DO
363# endif
364!
365!-----------------------------------------------------------------------
366! Compute tangent linear vertical mixing coefficients for momentum and
367! tracers. Compute S-coordinate vertical velocity, diagnostically from
368! horizontal mass divergence.
369!-----------------------------------------------------------------------
370!
371 DO ng=1,ngrids
372 DO tile=last_tile(ng),first_tile(ng),-1
373# if defined ANA_VMIX_NOT_YET
374 CALL rp_ana_vmix (ng, tile)
375# elif defined LMD_MIXING_NOT_YET
376 CALL rp_lmd_vmix (ng, tile)
377# elif defined BVF_MIXING_NOT_YET
378 CALL rp_bvf_mix (ng, tile)
379# endif
380 CALL rp_omega (ng, tile, irpm)
381!! CALL wvelocity (ng, tile, nstp(ng))
382 END DO
383!$OMP BARRIER
384 END DO
385!
386!-----------------------------------------------------------------------
387! Set free-surface to it time-averaged value. If applicable,
388! accumulate time-averaged output data which needs a irreversible
389! loop in shared-memory jobs.
390!-----------------------------------------------------------------------
391!
392 DO ng=1,ngrids
393 DO tile=first_tile(ng),last_tile(ng),+1 ! irreversible
394 CALL rp_set_zeta (ng, tile)
395# ifdef DIAGNOSTICS
396!! CALL rp_set_diags (ng, tile)
397# endif
398# ifdef RP_AVERAGES
399 CALL tl_set_avg (ng, tile)
400# endif
401 END DO
402!$OMP BARRIER
403 END DO
404!
405!-----------------------------------------------------------------------
406! If appropriate, write out fields into output NetCDF files. Notice
407! that IO data is written in delayed and serial mode. Exit if last
408! time step.
409!-----------------------------------------------------------------------
410!
411 DO ng=1,ngrids
412!$OMP MASTER
413 CALL rp_output (ng)
414!$OMP END MASTER
415!$OMP BARRIER
416 IF ((founderror(exit_flag, noerror, __line__, myfile)).or. &
417 & ((iic(ng).eq.(ntend(ng)+1)).and.(ng.eq.ngrids))) RETURN
418 END DO
419!
420!-----------------------------------------------------------------------
421! Compute right-hand-side terms for 3D equations.
422!-----------------------------------------------------------------------
423!
424 DO ng=1,ngrids
425 DO tile=first_tile(ng),last_tile(ng),+1
426 CALL rp_rhs3d (ng, tile)
427# ifdef MY25_MIXING_NOT_YET
428 CALL rp_my25_prestep (ng, tile)
429# elif defined GLS_MIXING_NOT_YET
430 CALL rp_gls_prestep (ng, tile)
431# endif
432 END DO
433!$OMP BARRIER
434 END DO
435!
436!-----------------------------------------------------------------------
437! Solve the vertically integrated primitive equations for the
438! free-surface and barotropic momentum components.
439!-----------------------------------------------------------------------
440!
441 loop_2d : DO my_iif=1,maxval(nfast)+1
442!
443! Set time indices for predictor step. The PREDICTOR_2D_STEP switch
444! it is assumed to be false before the first time-step.
445!
446 DO ng=1,ngrids
447 next_indx1=3-indx1(ng)
448 IF (.not.predictor_2d_step(ng).and. &
449 & my_iif.le.(nfast(ng)+1)) THEN
450 predictor_2d_step(ng)=.true.
451 iif(ng)=my_iif
452 IF (first_2d_step) THEN
453 kstp(ng)=indx1(ng)
454 ELSE
455 kstp(ng)=3-indx1(ng)
456 END IF
457 knew(ng)=3
458 krhs(ng)=indx1(ng)
459 END IF
460!
461! Predictor step - Advance barotropic equations using 2D time-step
462! ============== predictor scheme. No actual time-stepping is
463! performed during the auxiliary (nfast+1) time-step. It is needed
464! to finalize the fast-time averaging of 2D fields, if any, and
465! compute the new time-evolving depths.
466!
467 IF (my_iif.le.(nfast(ng)+1)) THEN
468 DO tile=last_tile(ng),first_tile(ng),-1
469 CALL rp_step2d (ng, tile)
470 END DO
471!$OMP BARRIER
472 END IF
473 END DO
474!
475! Set time indices for corrector step.
476!
477 DO ng=1,ngrids
478 IF (predictor_2d_step(ng)) THEN
479 predictor_2d_step(ng)=.false.
480 knew(ng)=next_indx1
481 kstp(ng)=3-knew(ng)
482 krhs(ng)=3
483 IF (iif(ng).lt.(nfast(ng)+1)) indx1(ng)=next_indx1
484 END IF
485!
486! Corrector step - Apply 2D time-step corrector scheme. Notice that
487! ============== there is not need for a corrector step during the
488! auxiliary (nfast+1) time-step.
489!
490 IF (iif(ng).lt.(nfast(ng)+1)) THEN
491 DO tile=first_tile(ng),last_tile(ng),+1
492 CALL rp_step2d (ng, tile)
493 END DO
494!$OMP BARRIER
495 END IF
496 END DO
497
498 END DO loop_2d
499!
500!-----------------------------------------------------------------------
501! Recompute depths and thicknesses using the new time filtered
502! free-surface.
503!-----------------------------------------------------------------------
504!
505 DO ng=1,ngrids
506 DO tile=first_tile(ng),last_tile(ng),+1
507 CALL rp_set_depth (ng, tile, irpm)
508 END DO
509!$OMP BARRIER
510 END DO
511!
512!-----------------------------------------------------------------------
513! Time-step 3D momentum equations.
514!-----------------------------------------------------------------------
515!
516! Time-step 3D momentum equations and couple with vertically
517! integrated equations.
518!
519 DO ng=1,ngrids
520 DO tile=last_tile(ng),first_tile(ng),-1
521 CALL rp_step3d_uv (ng, tile)
522 END DO
523!$OMP BARRIER
524 END DO
525!
526!-----------------------------------------------------------------------
527! Time-step vertical mixing turbulent equations and passive tracer
528! source and sink terms, if applicable.
529!-----------------------------------------------------------------------
530!
531 DO ng=1,ngrids
532 DO tile=first_tile(ng),last_tile(ng),+1
533 CALL rp_omega (ng, tile, irpm)
534# ifdef MY25_MIXING_NOT_YET
535 CALL rp_my25_corstep (ng, tile)
536# elif defined GLS_MIXING_NOT_YET
537 CALL rp_gls_corstep (ng, tile)
538# endif
539# ifdef BIOLOGY
540 CALL rp_biology (ng, tile)
541# endif
542# ifdef SEDIMENT_NOT_YET
543 CALL rp_sediment (ng, tile)
544# endif
545 END DO
546!$OMP BARRIER
547 END DO
548
549# ifndef TS_FIXED
550!
551!-----------------------------------------------------------------------
552! Time-step tracer equations.
553!-----------------------------------------------------------------------
554!
555 DO ng=1,ngrids
556 DO tile=last_tile(ng),first_tile(ng),-1
557 CALL rp_step3d_t (ng, tile)
558 END DO
559 END DO
560!$OMP BARRIER
561# endif
562
563# ifdef FLOATS_NOT_YET
564!
565!-----------------------------------------------------------------------
566! Compute Lagrangian drifters trajectories: Split all the drifters
567! between all the computational threads, except in distributed-memory
568! and serial configurations. In distributed-memory, the parallel node
569! containing the drifter is selected internally since the state
570! variables do not have a global scope.
571!-----------------------------------------------------------------------
572!
573 DO ng=1,ngrids
574 IF (lfloats(ng)) THEN
575# ifdef _OPENMP
576 chunk_size=(nfloats(ng)+numthreads-1)/numthreads
577 lstr=1+my_thread*chunk_size
578 lend=min(nfloats(ng),lstr+chunk_size-1)
579# else
580 lstr=1
581 lend=nfloats(ng)
582# endif
583 CALL rp_step_floats (ng, lstr, lend)
584!$OMP BARRIER
585!
586! Shift floats time indices.
587!
588 nfp1(ng)=mod(nfp1(ng)+1,nft+1)
589 nf(ng)=mod(nf(ng)+1,nft+1)
590 nfm1(ng)=mod(nfm1(ng)+1,nft+1)
591 nfm2(ng)=mod(nfm2(ng)+1,nft+1)
592 nfm3(ng)=mod(nfm3(ng)+1,nft+1)
593 END IF
594 END DO
595# endif
596 END DO step_loop
597
598 RETURN
599 END SUBROUTINE rp_main3d
600#else
601 SUBROUTINE rp_main3d
602 RETURN
603 END SUBROUTINE rp_main3d
604#endif
subroutine ana_vmix(ng, tile, model)
Definition ana_vmix.h:3
subroutine, public time_string(mytime, date_string)
Definition dateclock.F:1272
integer, dimension(:,:), allocatable couplesteps
integer iatmos
integer iwaves
integer numthreads
integer, dimension(:), allocatable first_tile
integer, dimension(:), allocatable last_tile
integer, dimension(:), allocatable nfloats
Definition mod_param.F:543
integer, parameter irpm
Definition mod_param.F:664
integer, parameter nft
Definition mod_param.F:539
integer ngrids
Definition mod_param.F:113
integer, parameter itlm
Definition mod_param.F:663
logical, dimension(:), allocatable lfloats
integer, dimension(:), allocatable iic
real(dp), dimension(:), allocatable dt
logical, dimension(:), allocatable predictor_2d_step
real(dp), dimension(:), allocatable tdays
logical, dimension(:), allocatable frequentimpulse
integer, dimension(:), allocatable nfast
real(dp), parameter sec2day
integer, dimension(:), allocatable ntend
character(len=22), dimension(:), allocatable time_code
integer exit_flag
integer, dimension(:), allocatable indx1
real(dp), dimension(:), allocatable time
integer, dimension(:), allocatable ntstart
integer nrun
integer, dimension(:), allocatable iif
integer noerror
integer, dimension(:), allocatable nfm2
integer, dimension(:), allocatable kstp
integer, dimension(:), allocatable knew
integer, dimension(:), allocatable nfm1
integer, dimension(:), allocatable nrhs
integer, dimension(:), allocatable lfinp
integer, dimension(:), allocatable lbinp
integer, dimension(:), allocatable nf
integer, dimension(:), allocatable nfm3
integer, dimension(:), allocatable nnew
integer, dimension(:), allocatable nfp1
integer, dimension(:), allocatable krhs
integer, dimension(:), allocatable nstp
subroutine, public omega(ng, tile, model)
Definition omega.F:42
subroutine, public rp_biology(ng, tile)
subroutine, public rp_bulk_flux(ng, tile)
subroutine, public rp_diag(ng, tile)
Definition rp_diag.F:26
subroutine, public rp_frc_adjust(ng, tile, linp)
subroutine, public rp_obc2d_adjust(ng, tile, linp)
subroutine, public rp_obc_adjust(ng, tile, linp)
subroutine, public rp_omega(ng, tile, model)
Definition rp_omega.F:35
subroutine, public rp_post_initial(ng, model)
subroutine, public rp_rho_eos(ng, tile, model)
Definition rp_rho_eos.F:48
subroutine, public rp_rhs3d(ng, tile)
Definition rp_rhs3d.F:29
subroutine, public rp_set_depth_bry(ng, tile, model)
subroutine, public rp_set_depth(ng, tile, model)
subroutine, public rp_set_massflux(ng, tile, model)
subroutine, public rp_set_vbc(ng, tile)
Definition rp_set_vbc.F:33
subroutine, public rp_set_zeta(ng, tile)
Definition rp_set_zeta.F:27
subroutine, public rp_step2d(ng, tile)
subroutine, public rp_step3d_t(ng, tile)
Definition rp_step3d_t.F:47
subroutine, public rp_step3d_uv(ng, tile)
subroutine, public set_depth(ng, tile, model)
Definition set_depth.F:34
subroutine, public set_massflux(ng, tile, model)
logical function, public founderror(flag, noerr, line, routine)
Definition strings.F:52
subroutine, public tl_forcing(ng, tile, kfrc, nfrc)
Definition tl_forcing.F:26
subroutine, public tl_set_avg(ng, tile)
Definition tl_set_avg.F:27
subroutine rp_get_data(ng)
Definition rp_get_data.F:4
subroutine rp_main3d(runinterval)
Definition rp_main3d.F:4
subroutine rp_output(ng)
Definition rp_output.F:4
subroutine rp_set_data(ng, tile)
Definition rp_set_data.F:4