source: flexpart.git/src/mpi_mod.f90 @ 54cbd6c

10.4.1_peseiGFS_025bugfixes+enhancementsdevrelease-10release-10.4.1scaling-bugunivie
Last change on this file since 54cbd6c was 79abee9, checked in by Espen Sollum ATMOS <eso@…>, 8 years ago

Fix for issue #151 (improper use of MPI_IN_PLACE). Also changed default RELEASES file to use SPECIES_001.

  • Property mode set to 100644
File size: 92.5 KB
Line 
1!**********************************************************************
2! Copyright 1998,1999,2000,2001,2002,2005,2007,2008,2009,2010         *
3! Andreas Stohl, Petra Seibert, A. Frank, Gerhard Wotawa,             *
4! Caroline Forster, Sabine Eckhardt, John Burkhart, Harald Sodemann   *
5!                                                                     *
6! This file is part of FLEXPART.                                      *
7!                                                                     *
8! FLEXPART is free software: you can redistribute it and/or modify    *
9! it under the terms of the GNU General Public License as published by*
10! the Free Software Foundation, either version 3 of the License, or   *
11! (at your option) any later version.                                 *
12!                                                                     *
13! FLEXPART is distributed in the hope that it will be useful,         *
14! but WITHOUT ANY WARRANTY; without even the implied warranty of      *
15! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the       *
16! GNU General Public License for more details.                        *
17!                                                                     *
18! You should have received a copy of the GNU General Public License   *
19! along with FLEXPART.  If not, see <http://www.gnu.org/licenses/>.   *
20!**********************************************************************
21
22module mpi_mod
23
24!*****************************************************************************
25!                                                                            *
26!  DESCRIPTION                                                               *
27!    This module contains subroutines and common variables used for the      *
28!    MPI parallelization of FLEXPART.                                        *
29!                                                                            *
30!  NOTE                                                                      *
31!    Depending on the MPI library installed on your system (e.g. mpich2,     *
32!    OpenMPI) you may need to choose below in this file between              *
33!      use mpi                                                               *
34!    (if the MPI library comes with the file 'mpi.mod'); or                  *
35!      include 'mpif.h'                                                      *
36!                                                                            *
37!                                                                            *
38!*****************************************************************************
39!                                                                            *
40! Variables:                                                                 *
41!                                                                            *
42! mp_ierr                 MPI error code                                     *
43! mp_np                   Number of MPI processes                            *
44! mp_pid                  Process ID of each MPI process                     *
45! mp_seed                 Parameter for random number seed                   *
46! read_grp_min            Minimum number of processes at which one will be   *
47!                         used as reader                                     *
48! numpart_mpi,            Number of particles per node                       *
49! maxpart_mpi                                                                *
50! mp_partid               MPI process ID for particle calculation            *
51! mp_partgroup_           Refers to the subset of processors performing      *
52!                         loops over particles. Will be all processes        *
53!                         unless a dedicated process runs getfields/readwind *
54! lmp_sync                If .false., use asynchronous MPI                   *
55! mp_cp                   Real precision to use for deposition fields        *
56!                                                                            *
57!                                                                            *
58!                                                                            *
59!                                                                            *
60!*****************************************************************************
61
62  use mpi
63  use par_mod, only: dp,sp
64  use com_mod, only: lroot
65
66  implicit none
67
68!  include 'mpif.h'
69
70  public
71
72! Set aside a process for reading windfields if using at least these many processes
73!==================================================
74  integer, parameter, private :: read_grp_min=4
75!==================================================
76
77! Variables for each MPI process in the world group
78  integer :: mp_ierr, mp_np, mp_pid, mp_partid
79  integer, private :: world_group_id
80
81! Variables for MPI processes in the 'particle' group
82  integer, allocatable, dimension(:) :: mp_partgroup_rank
83  integer :: mp_partgroup_comm, mp_partgroup_pid, mp_partgroup_np
84
85  integer :: mp_seed=0
86  integer, parameter :: mp_sp=MPI_REAL4, mp_dp=MPI_REAL8
87  integer :: mp_cp
88  integer, parameter :: id_root=0 ! master process
89
90! MPI tags/requests for send/receive operation
91  integer :: tm1
92  integer, parameter :: nvar_async=26 !27 !29 :DBG:
93!integer, dimension(:), allocatable :: tags
94  integer, dimension(:), allocatable :: reqs
95
96
97  integer :: id_read   ! readwind/getfield process
98  integer :: numpart_mpi,maxpart_mpi ! number of particles per node
99  integer :: tot_numpart=0
100  integer :: mp_comm_used ! global or subgroup communicator
101
102  logical :: lmpreader=.false. ! is set to true for reading process(es) only.
103  logical :: lmp_use_reader=.false. ! true if separate readwind process is used
104
105! .true. if only using synchronous MPI send/recv (default)
106! If setting this to .false., numwfmem must be set to 3
107!===============================================================================
108  logical :: lmp_sync=.true.
109!===============================================================================
110
111! mp_dbg_mode       Used for debugging MPI.
112! mp_dev_mode       various checks related to debugging the parallel code
113! mp_dbg_out        write some arrays to data file for debugging
114! mp_measure_time   Measure cpu/wall time, write out at end of run
115! mp_time_barrier   Measure MPI barrier time
116! mp_exact_numpart  Use an extra MPI communication to give the exact number of particles
117!                   to standard output (this does *not* otherwise affect the simulation)
118  logical, parameter :: mp_dbg_mode = .false.
119  logical, parameter :: mp_dev_mode = .false.
120  logical, parameter :: mp_dbg_out = .false.
121  logical, parameter :: mp_time_barrier=.true.
122  logical, parameter :: mp_measure_time=.false.
123  logical, parameter :: mp_exact_numpart=.true.
124
125! for measuring CPU/Wall time
126  real(sp),private :: mp_comm_time_beg, mp_comm_time_end, mp_comm_time_total=0.
127  real(dp),private :: mp_comm_wtime_beg, mp_comm_wtime_end, mp_comm_wtime_total=0.
128  real(sp),private :: mp_root_time_beg, mp_root_time_end, mp_root_time_total=0.
129  real(dp),private :: mp_root_wtime_beg, mp_root_wtime_end, mp_root_wtime_total=0.
130  real(sp),private :: mp_barrier_time_beg, mp_barrier_time_end, mp_barrier_time_total=0.
131  real(dp),private :: mp_barrier_wtime_beg, mp_barrier_wtime_end, mp_barrier_wtime_total=0.
132  real(sp),private :: tm_nploop_beg, tm_nploop_end, tm_nploop_total=0.
133  real(sp),private :: tm_tot_beg, tm_tot_end, tm_tot_total=0.
134  real(dp),private :: mp_getfields_wtime_beg, mp_getfields_wtime_end, mp_getfields_wtime_total=0.
135  real(sp),private :: mp_getfields_time_beg, mp_getfields_time_end, mp_getfields_time_total=0.
136  real(dp),private :: mp_readwind_wtime_beg, mp_readwind_wtime_end, mp_readwind_wtime_total=0.
137  real(sp),private :: mp_readwind_time_beg, mp_readwind_time_end, mp_readwind_time_total=0.
138  real(dp),private :: mp_io_wtime_beg, mp_io_wtime_end, mp_io_wtime_total=0.
139  real(sp),private :: mp_io_time_beg, mp_io_time_end, mp_io_time_total=0.
140  real(dp),private :: mp_wetdepo_wtime_beg, mp_wetdepo_wtime_end, mp_wetdepo_wtime_total=0.
141  real(sp),private :: mp_wetdepo_time_beg, mp_wetdepo_time_end, mp_wetdepo_time_total=0.
142  real(dp),private :: mp_advance_wtime_beg, mp_advance_wtime_end, mp_advance_wtime_total=0.
143  real(dp),private :: mp_conccalc_time_beg, mp_conccalc_time_end, mp_conccalc_time_total=0.
144  real(dp),private :: mp_total_wtime_beg, mp_total_wtime_end, mp_total_wtime_total=0.
145  real(dp),private :: mp_vt_wtime_beg, mp_vt_wtime_end, mp_vt_wtime_total
146  real(sp),private :: mp_vt_time_beg, mp_vt_time_end, mp_vt_time_total
147
148! dat_lun           logical unit number for i/o
149  integer, private :: dat_lun
150
151contains
152
153  subroutine mpif_init
154!***********************************************************************
155! DESCRIPTION
156!   Initialize MPI.
157!   
158!   Create the global communicator MPI_COMM_WORLD
159!   If using a separate MPI process for getfields/readwind, a subgroup
160!   is created for the other processes.
161!
162! VARIABLES
163!   mpi_mode    default 0, set to 2/3 if running MPI version
164!   mp_np       number of running processes, decided at run-time
165!***********************************************************************
166    use par_mod, only: maxpart, numwfmem, dep_prec
167    use com_mod, only: mpi_mode, verbosity
168
169    implicit none
170
171    integer :: i,j,s,addmaxpart=0
172
173! Each process gets an ID (mp_pid) in the range 0,..,mp_np-1
174    call MPI_INIT(mp_ierr)
175    if (mp_ierr /= 0) goto 100
176    call MPI_COMM_RANK(MPI_COMM_WORLD, mp_pid, mp_ierr)
177    if (mp_ierr /= 0) goto 100
178    call MPI_COMM_SIZE(MPI_COMM_WORLD, mp_np, mp_ierr)
179    if (mp_ierr /= 0) goto 100
180
181
182! Variable mpi_mode is used to handle subroutines common to parallel/serial version
183    if (lmp_sync) then
184      mpi_mode=2 ! hold 2 windfields in memory
185    else
186      mpi_mode=3 ! hold 3 windfields in memory
187    end if
188
189    if (mp_pid.ne.0) then
190      lroot = .false.
191    end if
192
193! Set MPI precision to use for transferring deposition fields
194!************************************************************
195    if (dep_prec==dp) then
196      mp_cp = MPI_REAL8
197      ! TODO: write info message for serial version as well
198      if (lroot.and.verbosity>0) write(*,*) 'Using double precision for deposition fields'
199    else if (dep_prec==sp) then
200      mp_cp = MPI_REAL4
201      if (lroot.and.verbosity>0) write(*,*) 'Using single precision for deposition fields'
202    else
203      write(*,*) 'ERROR: something went wrong setting MPI real precision'
204      stop
205    end if
206
207! Check for sensible combination of parameters
208!*********************************************
209
210    if (.not.lmp_sync.and.numwfmem.ne.3) then
211      if (lroot) then
212        write(*,FMT='(80("#"))')
213        write(*,*) '#### mpi_mod::mpif_init> ERROR: ', &
214             & 'numwfmem must be set to 3 for asyncronous reading ####'
215        write(*,FMT='(80("#"))')
216      end if
217      call MPI_FINALIZE(mp_ierr)
218      stop
219    else if (lmp_sync.and.numwfmem.ne.2.and.lroot) then
220      write(*,FMT='(80("#"))')
221      write(*,*) '#### mpi_mod::mpif_init> WARNING: ', &
222           & 'numwfmem should be set to 2 for syncronous'
223      write(*,*) ' reading. Results will still be valid, but unneccesary memory &
224           &is allocated.'
225      write(*,FMT='(80("#"))')
226! Force "syncronized" version if all processes will call getfields
227    else if ((.not.lmp_sync.and.mp_np.lt.read_grp_min).or.(mp_np.eq.1)) then
228      if (lroot) then
229        write(*,FMT='(80("#"))')
230        write(*,*) '#### mpi_mod::mpif_init> WARNING: ', &
231             & 'all procs call getfields. Setting lmp_sync=.true.'
232        write(*,FMT='(80("#"))')
233      end if
234      lmp_sync=.true. ! :DBG: eso fix this...
235    end if
236
237! TODO: Add more warnings for unimplemented flexpart features
238
239! Set ID of process that calls getfield/readwind.
240! Using the last process in the group ensures statistical identical results
241! as running with one process less but not using separate read process
242!**********************************************************************
243
244!    id_read = min(mp_np-1, 1)
245    id_read = mp_np-1
246
247    if (mp_pid.eq.id_read) lmpreader=.true.
248
249    call MPI_Comm_group (MPI_COMM_WORLD, world_group_id, mp_ierr)
250
251! Create a MPI group/communicator that will calculate trajectories.
252! Skip this step if program is run with only a few processes
253!************************************************************************
254    allocate(mp_partgroup_rank(0:mp_np-2))
255
256! This allows checking for allocation of arrays when no subgroup is used
257    mp_partgroup_pid=0
258
259    if (read_grp_min.lt.2) then
260      write(*,*) '#### mpi_mod::mpif_init> ERROR ####', &
261           & 'read_grp_min should be at least 2. Exiting'
262      stop
263    end if
264
265    if (mp_np.ge.read_grp_min) then
266      lmp_use_reader = .true.
267
268! Map the subgroup IDs= 0,1,2,...,mp_np-2, skipping reader process
269      j=-1
270      do i=0, mp_np-2 !loop over all (subgroup) IDs
271        j=j+1
272        if (i.eq.id_read) then
273          j=j+1
274        end if
275        mp_partgroup_rank(i) = j
276      end do
277
278      call MPI_Group_incl (world_group_id, mp_np-1, mp_partgroup_rank, &
279           &mp_partgroup_pid, mp_ierr)
280      if (mp_ierr /= 0) goto 100
281      call MPI_Comm_create (MPI_COMM_WORLD, mp_partgroup_pid, mp_partgroup_comm, mp_ierr)
282      if (mp_ierr /= 0) goto 100
283
284      if (mp_pid.ne.id_read) then
285        call MPI_Comm_rank (mp_partgroup_comm, mp_partgroup_pid, mp_ierr)
286        if (mp_ierr /= 0) goto 100
287
288! Get size of new group (=mp_np-1)
289        call MPI_COMM_SIZE(mp_partgroup_comm, mp_partgroup_np, mp_ierr)
290        if (mp_ierr /= 0) goto 100
291        if (mp_partgroup_np.ne.mp_np-1) then
292          write(*,*)  '#### mpi_mod:: mpif_init> ERROR ####&
293               & mp_partgroup_np.ne.mp_np-1'
294          stop
295        endif
296
297      else
298        mp_partgroup_pid = -1
299      end if
300    end if
301
302    if (lmp_use_reader) then
303      mp_comm_used = mp_partgroup_comm
304      mp_partid = mp_partgroup_pid
305      mp_partgroup_np=mp_np-1
306    else
307      mp_comm_used = MPI_COMM_WORLD
308      mp_partgroup_np = mp_np
309      mp_partid = mp_pid
310    end if
311
312! Set maxpart per process
313    if (mp_partid.lt.mod(maxpart,mp_partgroup_np)) addmaxpart=1
314    maxpart_mpi=int(maxpart/mp_partgroup_np)+addmaxpart
315
316! Do not allocate particle data arrays for readwind process
317    if (lmpreader.and.lmp_use_reader) then
318      maxpart_mpi=0
319    end if
320
321! Set random seed for each non-root process
322    if (mp_pid.gt.0) then
323!    if (mp_pid.ge.0) then
324!      call system_clock(s)
325      s = 244
326      mp_seed = -abs(mod((s*181)*((mp_pid-83)*359), 104729))
327    end if
328    if (mp_dev_mode) write(*,*) 'PID, mp_seed: ',mp_pid, mp_seed
329    if (mp_dbg_mode) then
330! :DBG: For debugging, set all seed to 0 and maxrand to e.g. 20
331      mp_seed=0
332      if (lroot) write(*,*) 'MPI: setting seed=0'
333    end if
334
335! Allocate request array for asynchronous MPI
336    if (.not.lmp_sync) then
337      allocate(reqs(0:nvar_async*mp_np-1))
338      reqs(:)=MPI_REQUEST_NULL
339    else
340      allocate(reqs(0:1))
341      reqs(:)=MPI_REQUEST_NULL
342    end if
343
344    goto 101
345
346100 write(*,*) '#### mpi_mod::mpif_init> ERROR ####', mp_ierr
347    stop
348
349101 end subroutine mpif_init
350
351
352  subroutine mpif_mpi_barrier
353!***********************************************************************
354! Set MPI barrier (all processes are synchronized here), with the
355! possible exception of a separate 'readwind' process.
356! Optionally measure cpu/wall time.
357!
358!***********************************************************************
359    implicit none
360
361    if (mp_time_barrier) then
362      call cpu_time(mp_barrier_time_beg)
363      mp_barrier_wtime_beg = mpi_wtime()
364    endif
365
366    call MPI_BARRIER(mp_comm_used, mp_ierr)
367    if (mp_ierr /= 0) goto 100
368
369    if (mp_time_barrier) then
370      call cpu_time(mp_barrier_time_end)
371      mp_barrier_wtime_end = mpi_wtime()
372
373      mp_barrier_time_total=mp_barrier_time_total+&
374           &(mp_barrier_time_end-mp_barrier_time_beg)
375      mp_barrier_wtime_total=mp_barrier_wtime_total+&
376           &(mp_barrier_wtime_end-mp_barrier_wtime_beg)
377    end if
378
379    goto 101
380
381100 write(*,*) '#### mpi_mod::mpif_mpi_barrier> ERROR ####', mp_ierr
382    stop
383
384101 end subroutine mpif_mpi_barrier
385
386
387  subroutine mpif_com_mod_allocate
388!*******************************************************************************   
389! Dynamic allocation of arrays (in serial code these have size maxpart)
390!
391! Not in use anymore, moved to com_mod for interoperability with serial version
392!
393! TODO: remove
394!*******************************************************************************
395    use com_mod
396
397    implicit none
398
399    integer :: array_size
400
401    array_size = maxpart_mpi
402
403    allocate(itra1(array_size),npoint(array_size),&
404         & nclass(array_size),idt(array_size),&
405         & itramem(array_size),itrasplit(array_size),&
406         & xtra1(array_size),ytra1(array_size),&
407         & ztra1(array_size),xmass1(array_size, maxspec))
408
409    allocate(uap(array_size),ucp(array_size),&
410         & uzp(array_size),us(array_size),&
411         & vs(array_size),&
412         & ws(array_size),cbt(array_size))
413
414
415  end subroutine mpif_com_mod_allocate
416
417
418  subroutine mpif_tm_send_dims
419!***********************************************************************
420! Distribute array dimensions from pid0 to all processes.
421! numpart_mpi/numpart is sent to allow dynamic allocation
422!
423! Currently not in use.
424!
425! Variables of similar type (integer, real) are placed in an array
426! and sent collectively, to avoid the overhead associated with individual
427! MPI send/receives
428!
429!
430!***********************************************************************
431    use com_mod
432
433    implicit none
434
435    integer :: add_part=0
436
437    call MPI_Bcast(numpart, 1, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
438
439! MPI subgroup does the particle-loop
440    if (lmp_use_reader) then
441      if (mod(numpart,mp_partgroup_np).ne.0) add_part=1
442      numpart_mpi=int(numpart/mp_partgroup_np)+add_part
443    else
444      if (mod(numpart,mp_np).ne.0) add_part=1
445      numpart_mpi=int(numpart/mp_np)+add_part
446    end if
447
448
449! redefine numpart as 'numpart per process' throughout the code
450!**************************************************************
451    numpart = numpart_mpi
452
453  end subroutine mpif_tm_send_dims
454
455
456  subroutine mpif_tm_send_vars
457!***********************************************************************
458! Distribute particle variables from pid0 to all processes.
459! Called from timemanager
460! *NOT IN USE* at the moment, but can be useful for debugging
461!
462!***********************************************************************
463    use com_mod
464
465    implicit none
466
467    integer :: i
468
469! Time for MPI communications
470!****************************
471    if (mp_measure_time) call mpif_mtime('commtime',0)
472
473! Distribute variables, send from pid 0 to other processes
474!**********************************************************************
475
476! integers
477    if (lroot) then
478      call MPI_SCATTER(npoint,numpart_mpi,MPI_INTEGER,MPI_IN_PLACE,&
479           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
480      if (mp_ierr /= 0) goto 600
481      call MPI_SCATTER(idt,numpart_mpi,MPI_INTEGER,MPI_IN_PLACE,&
482           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
483      if (mp_ierr /= 0) goto 600
484      call MPI_SCATTER(itra1,numpart_mpi,MPI_INTEGER,MPI_IN_PLACE,&
485           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
486      if (mp_ierr /= 0) goto 600
487      call MPI_SCATTER(nclass,numpart_mpi,MPI_INTEGER,MPI_IN_PLACE,&
488           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
489      if (mp_ierr /= 0) goto 600
490      call MPI_SCATTER(itramem,numpart_mpi,MPI_INTEGER,MPI_IN_PLACE,&
491           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
492      if (mp_ierr /= 0) goto 600
493
494! int2
495      call MPI_SCATTER(cbt,numpart_mpi,MPI_INTEGER2,MPI_IN_PLACE,&
496           &numpart_mpi,MPI_INTEGER2,id_root,mp_comm_used,mp_ierr)
497      if (mp_ierr /= 0) goto 600
498
499! real
500      call MPI_SCATTER(uap,numpart_mpi,mp_sp,MPI_IN_PLACE,&
501           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
502      if (mp_ierr /= 0) goto 600
503      call MPI_SCATTER(ucp,numpart_mpi,mp_sp,MPI_IN_PLACE,&
504           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
505      if (mp_ierr /= 0) goto 600
506      call MPI_SCATTER(uzp,numpart_mpi,mp_sp,MPI_IN_PLACE,&
507           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
508      if (mp_ierr /= 0) goto 600
509
510      call MPI_SCATTER(us,numpart_mpi,mp_sp,MPI_IN_PLACE,&
511           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
512      if (mp_ierr /= 0) goto 600
513      call MPI_SCATTER(vs,numpart_mpi,mp_sp,MPI_IN_PLACE,&
514           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
515      if (mp_ierr /= 0) goto 600
516      call MPI_SCATTER(ws,numpart_mpi,mp_sp,MPI_IN_PLACE,&
517           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
518      if (mp_ierr /= 0) goto 600
519
520      call MPI_SCATTER(xtra1,numpart_mpi,mp_dp,MPI_IN_PLACE,&
521           &numpart_mpi,mp_dp,id_root,mp_comm_used,mp_ierr)
522      if (mp_ierr /= 0) goto 600
523      call MPI_SCATTER(ytra1,numpart_mpi,mp_dp,MPI_IN_PLACE,&
524           &numpart_mpi,mp_dp,id_root,mp_comm_used,mp_ierr)
525      if (mp_ierr /= 0) goto 600
526      call MPI_SCATTER(ztra1,numpart_mpi,mp_sp,MPI_IN_PLACE,&
527           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
528      if (mp_ierr /= 0) goto 600
529
530      do i=1,nspec
531        call MPI_SCATTER(xmass1(:,i),numpart_mpi,mp_sp,MPI_IN_PLACE,&
532             &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
533        if (mp_ierr /= 0) goto 600
534      end do
535
536    else ! (mp_pid >= 1)
537! integers
538      call MPI_SCATTER(npoint,numpart_mpi,MPI_INTEGER,npoint,&
539           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
540      if (mp_ierr /= 0) goto 600
541      call MPI_SCATTER(idt,numpart_mpi,MPI_INTEGER,idt,&
542           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
543      if (mp_ierr /= 0) goto 600
544      call MPI_SCATTER(itra1,numpart_mpi,MPI_INTEGER,itra1,&
545           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
546      if (mp_ierr /= 0) goto 600
547      call MPI_SCATTER(nclass,numpart_mpi,MPI_INTEGER,nclass,&
548           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
549      if (mp_ierr /= 0) goto 600
550      call MPI_SCATTER(itramem,numpart_mpi,MPI_INTEGER,itramem,&
551           &numpart_mpi,MPI_INTEGER,id_root,mp_comm_used,mp_ierr)
552      if (mp_ierr /= 0) goto 600
553
554! int2
555      call MPI_SCATTER(cbt,numpart_mpi,MPI_INTEGER2,cbt,&
556           &numpart_mpi,MPI_INTEGER2,id_root,mp_comm_used,mp_ierr)
557      if (mp_ierr /= 0) goto 600
558
559! reals
560      call MPI_SCATTER(uap,numpart_mpi,mp_sp,uap,&
561           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
562      if (mp_ierr /= 0) goto 600
563      call MPI_SCATTER(ucp,numpart_mpi,mp_sp,ucp,&
564           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
565      if (mp_ierr /= 0) goto 600
566      call MPI_SCATTER(uzp,numpart_mpi,mp_sp,uzp,&
567           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
568      if (mp_ierr /= 0) goto 600
569
570      call MPI_SCATTER(us,numpart_mpi,mp_sp,us,&
571           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
572      if (mp_ierr /= 0) goto 600
573      call MPI_SCATTER(vs,numpart_mpi,mp_sp,vs,&
574           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
575      if (mp_ierr /= 0) goto 600
576      call MPI_SCATTER(ws,numpart_mpi,mp_sp,ws,&
577           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
578      if (mp_ierr /= 0) goto 600
579
580      call MPI_SCATTER(xtra1,numpart_mpi,mp_dp,xtra1,&
581           &numpart_mpi,mp_dp,id_root,mp_comm_used,mp_ierr)
582      if (mp_ierr /= 0) goto 600
583      call MPI_SCATTER(ytra1,numpart_mpi,mp_dp,ytra1,&
584           &numpart_mpi,mp_dp,id_root,mp_comm_used,mp_ierr)
585      if (mp_ierr /= 0) goto 600
586      call MPI_SCATTER(ztra1,numpart_mpi,mp_sp,ztra1,&
587           &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
588      if (mp_ierr /= 0) goto 600
589
590      do i=1,nspec
591        call MPI_SCATTER(xmass1(:,i),numpart_mpi,mp_sp,xmass1(:,i),&
592             &numpart_mpi,mp_sp,id_root,mp_comm_used,mp_ierr)
593        if (mp_ierr /= 0) goto 600
594      end do
595
596    end if !(lroot)
597
598    if (mp_measure_time) call mpif_mtime('commtime',1)
599
600    goto 601
601
602600 write(*,*) "mpi_mod> mp_ierr \= 0", mp_ierr
603    stop
604601 end subroutine mpif_tm_send_vars
605
606
607  subroutine mpif_tm_collect_vars
608!*******************************************************************************
609! Collect particle data to PID 0 from all processes.                           *
610!                                                                              *
611! This can be called from timemanager_mpi, before caclulating                  *
612! concentrations/writing output grids.                                         *
613!                                                                              *
614! Currently *not in use*, as each process calculates concentrations            *
615! separately, but can be useful for debugging                                  *
616!                                                                              *
617! To use this routine (together with mpif_tm_send_vars) to transfer data       *
618! to the MPI root process (useful for debugging), insert code like this        *
619! (in timemanager_mpi.f90)                                                     *
620!                                                                              *
621!      if (lroot) tot_numpart=numpart ! root stores total numpart              *
622!      call mpif_tm_send_dims                                                  *
623!      if (numpart>1) then                                                     *
624!        call mpif_tm_send_vars                                                *
625!      end if                                                                  *
626!                                                                              *
627!    To collect data from all processes to root process after a parallel       *
628!    region:                                                                   *
629!                                                                              *
630!      if (numpart>0) then                                                     *
631!          if (lroot) numpart = tot_numpart                                    *
632!       call mpif_tm_collect_vars                                              *
633!      end if                                                                  *
634!                                                                              *
635! Then a section that begins with "if (lroot) ..." will behave like            *
636! serial flexpart                                                              *
637!                                                                              *
638!*******************************************************************************
639    use com_mod !, only: numpart, nspec, itra1, npoint, nclass
640
641    implicit none
642
643    integer :: i
644
645    logical :: use_mp_vars = .false.! use mp_... type allocated vars
646    logical :: use_in_place = .true.! use MPI_IN_PLACE to collect vars.
647
648
649! Time for MPI communications
650    if (mp_measure_time) call mpif_mtime('commtime',0)
651
652! Distribute variables, send from pid 0 to other processes
653!**********************************************************************
654
655    if (.not. use_mp_vars.and..not.use_in_place) then
656
657! Integers:
658      call MPI_GATHER(npoint, numpart_mpi, MPI_INTEGER, npoint, &
659           &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
660      if (mp_ierr /= 0) goto 600
661      call MPI_GATHER(idt, numpart_mpi, MPI_INTEGER, idt, &
662           &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
663      if (mp_ierr /= 0) goto 600
664      call MPI_GATHER(itra1, numpart_mpi, MPI_INTEGER, itra1, &
665           &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
666      if (mp_ierr /= 0) goto 600
667      call MPI_GATHER(nclass, numpart_mpi, MPI_INTEGER, nclass, &
668           &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
669      if (mp_ierr /= 0) goto 600
670      call MPI_GATHER(itramem, numpart_mpi, MPI_INTEGER, itramem, &
671           &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
672      if (mp_ierr /= 0) goto 600
673
674! int2
675      call MPI_GATHER(cbt, numpart_mpi, MPI_INTEGER2, cbt, &
676           &numpart_mpi, MPI_INTEGER2, id_root, mp_comm_used, mp_ierr)
677      if (mp_ierr /= 0) goto 600
678
679! Reals:
680      call MPI_GATHER(uap, numpart_mpi, mp_sp, uap, &
681           &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
682      if (mp_ierr /= 0) goto 600
683      call MPI_GATHER(ucp, numpart_mpi, mp_sp, ucp, &
684           &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
685      if (mp_ierr /= 0) goto 600
686      call MPI_GATHER(uzp, numpart_mpi, mp_sp, uzp, &
687           &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
688      if (mp_ierr /= 0) goto 600
689
690      call MPI_GATHER(us, numpart_mpi, mp_sp, us, &
691           &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
692      if (mp_ierr /= 0) goto 600
693      call MPI_GATHER(vs, numpart_mpi, mp_sp, vs, &
694           &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
695      if (mp_ierr /= 0) goto 600
696      call MPI_GATHER(ws, numpart_mpi, mp_sp, ws, &
697           &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
698      if (mp_ierr /= 0) goto 600
699
700
701      call MPI_GATHER(xtra1, numpart_mpi, mp_dp, xtra1, &
702           &numpart_mpi, mp_dp, id_root, mp_comm_used, mp_ierr)
703      if (mp_ierr /= 0) goto 600
704      call MPI_GATHER(ytra1, numpart_mpi, mp_dp, ytra1, &
705           &numpart_mpi, mp_dp, id_root, mp_comm_used, mp_ierr)
706      if (mp_ierr /= 0) goto 600
707      call MPI_GATHER(ztra1, numpart_mpi, mp_sp, ztra1, &
708           &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
709      if (mp_ierr /= 0) goto 600
710
711      do i=1, nspec
712        call MPI_GATHER(xmass1(:,i), numpart_mpi, mp_sp, xmass1(:,i), &
713             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
714        if (mp_ierr /= 0) goto 600
715
716      end do
717
718! Use variables npoint etc directly for communications
719!***********************************************************************
720    else if (use_in_place.and..not.use_mp_vars) then
721      if (lroot) then
722! Integers:
723        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, MPI_INTEGER, npoint, &
724             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
725        if (mp_ierr /= 0) goto 600
726        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, MPI_INTEGER, idt, &
727             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
728        if (mp_ierr /= 0) goto 600
729        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, MPI_INTEGER, itra1, &
730             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
731        if (mp_ierr /= 0) goto 600
732        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, MPI_INTEGER, nclass, &
733             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
734        if (mp_ierr /= 0) goto 600
735        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, MPI_INTEGER, itramem, &
736             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
737        if (mp_ierr /= 0) goto 600
738
739! int2
740        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, MPI_INTEGER2, cbt, &
741             &numpart_mpi, MPI_INTEGER2, id_root, mp_comm_used, mp_ierr)
742        if (mp_ierr /= 0) goto 600
743
744! Reals:
745        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, uap, &
746             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
747        if (mp_ierr /= 0) goto 600
748        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, ucp, &
749             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
750        if (mp_ierr /= 0) goto 600
751        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, uzp, &
752             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
753        if (mp_ierr /= 0) goto 600
754
755        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, us, &
756             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
757        if (mp_ierr /= 0) goto 600
758        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, vs, &
759             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
760        if (mp_ierr /= 0) goto 600
761        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, ws, &
762             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
763        if (mp_ierr /= 0) goto 600
764
765
766        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_dp, xtra1, &
767             &numpart_mpi, mp_dp, id_root, mp_comm_used, mp_ierr)
768        if (mp_ierr /= 0) goto 600
769        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_dp, ytra1, &
770             &numpart_mpi, mp_dp, id_root, mp_comm_used, mp_ierr)
771        if (mp_ierr /= 0) goto 600
772        call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, ztra1, &
773             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
774        if (mp_ierr /= 0) goto 600
775
776        do i=1, nspec
777          call MPI_GATHER(MPI_IN_PLACE, numpart_mpi, mp_sp, xmass1(:,i), &
778               &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
779          if (mp_ierr /= 0) goto 600
780        end do
781
782      else ! (for mp_pid >= 1)
783
784! Integers:
785        call MPI_GATHER(npoint, numpart_mpi, MPI_INTEGER, npoint, &
786             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
787        if (mp_ierr /= 0) goto 600
788        call MPI_GATHER(idt, numpart_mpi, MPI_INTEGER, idt, &
789             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
790        if (mp_ierr /= 0) goto 600
791        call MPI_GATHER(itra1, numpart_mpi, MPI_INTEGER, itra1, &
792             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
793        if (mp_ierr /= 0) goto 600
794        call MPI_GATHER(nclass, numpart_mpi, MPI_INTEGER, nclass, &
795             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
796        if (mp_ierr /= 0) goto 600
797        call MPI_GATHER(itramem, numpart_mpi, MPI_INTEGER, itramem, &
798             &numpart_mpi, MPI_INTEGER, id_root, mp_comm_used, mp_ierr)
799        if (mp_ierr /= 0) goto 600
800
801! int2
802        call MPI_GATHER(cbt, numpart_mpi, MPI_INTEGER2, cbt, &
803             &numpart_mpi, MPI_INTEGER2, id_root, mp_comm_used, mp_ierr)
804        if (mp_ierr /= 0) goto 600
805
806! Reals:
807        call MPI_GATHER(uap, numpart_mpi, mp_sp, uap, &
808             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
809        if (mp_ierr /= 0) goto 600
810        call MPI_GATHER(ucp, numpart_mpi, mp_sp, ucp, &
811             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
812        if (mp_ierr /= 0) goto 600
813        call MPI_GATHER(uzp, numpart_mpi, mp_sp, uzp, &
814             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
815        if (mp_ierr /= 0) goto 600
816
817        call MPI_GATHER(us, numpart_mpi, mp_sp, us, &
818             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
819        if (mp_ierr /= 0) goto 600
820        call MPI_GATHER(vs, numpart_mpi, mp_sp, vs, &
821             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
822        if (mp_ierr /= 0) goto 600
823        call MPI_GATHER(ws, numpart_mpi, mp_sp, ws, &
824             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
825        if (mp_ierr /= 0) goto 600
826
827
828        call MPI_GATHER(xtra1, numpart_mpi, mp_dp, xtra1, &
829             &numpart_mpi, mp_dp, id_root, mp_comm_used, mp_ierr)
830        if (mp_ierr /= 0) goto 600
831        call MPI_GATHER(ytra1, numpart_mpi, mp_dp, ytra1, &
832             &numpart_mpi, mp_dp, id_root, mp_comm_used, mp_ierr)
833        if (mp_ierr /= 0) goto 600
834        call MPI_GATHER(ztra1, numpart_mpi, mp_sp, ztra1, &
835             &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
836        if (mp_ierr /= 0) goto 600
837
838        do i=1, nspec
839          call MPI_GATHER(xmass1(:,i), numpart_mpi, mp_sp, xmass1(:,i), &
840               &numpart_mpi, mp_sp, id_root, mp_comm_used, mp_ierr)
841          if (mp_ierr /= 0) goto 600
842        end do
843      end if ! (lroot)
844    end if !  (use_in_place)
845
846    if (mp_measure_time) call mpif_mtime('commtime',1)
847
848    goto 601
849
850600 write(*,*) "mpi_mod> mp_ierr \= 0", mp_ierr
851    stop
852601 end subroutine mpif_tm_collect_vars
853
854
855  subroutine mpif_gf_send_vars(memstat)
856!*******************************************************************************
857! DESCRIPTION
858!   Distribute 'getfield' variables from reader process
859!
860!   Called from timemanager
861!
862! NOTE
863!   This subroutine distributes windfields read from the reader process to
864!   all other processes. Usually one set of fields are transfered, but at
865!   first timestep when there are no fields in memory, both are transfered.
866!   MPI_Bcast is used, so implicitly all processes are synchronized at this
867!   step
868!
869!
870!*******************************************************************************
871    use com_mod
872    use par_mod,only: numwfmem
873
874    implicit none
875
876    integer, intent(in) :: memstat
877
878! Common array sizes used for communications
879    integer, parameter :: d3_size1 = nxmax*nymax*nzmax
880    integer, parameter :: d3_size2 = nxmax*nymax*nuvzmax
881    integer, parameter :: d2_size1 = nxmax*nymax
882    integer, parameter :: d2_size2 = nxmax*nymax*maxspec
883    integer, parameter :: d2_size3 = nxmax*nymax
884    integer, parameter :: d1_size1 = maxwf
885
886    integer :: d3s1,d3s2,d2s1,d2s2
887
888! li,ui    lower,upper indices of fields to be transfered (1-3, ui>=li)
889    integer :: li,ui
890
891! First time routine is called the unchangeable fields will be transfered   
892    logical :: first_call=.true.
893
894!**********************************************************************
895
896! Sizes of arrays transfered
897    d3s1=d3_size1
898    d3s2=d3_size2
899    d2s1=d2_size1
900    d2s2=d2_size2
901
902! Decide which field(s) need to be transfered
903    if (memstat.ge.32) then ! distribute 2 fields, to 1st/2nd indices
904      li=1; ui=2
905      d3s1=2*d3_size1
906      d3s2=2*d3_size2
907      d2s1=2*d2_size1
908      d2s2=2*d2_size2
909    else if (memstat.eq.17) then ! distribute 1 field, place on 1st index
910      li=1; ui=1
911    else if (memstat.eq.18) then ! distribute 1 field, place on 2nd index
912      li=2; ui=2
913    else if (memstat.eq.19) then ! distribute 1 field, place on 3nd index
914      li=3; ui=3
915    else
916      write(*,*) '#### mpi_mod::mpif_gf_send_vars> ERROR: ', &
917           & 'wrong value of memstat, exiting ####', memstat
918      stop
919    end if
920
921
922! Time for MPI communications
923    if (mp_measure_time) call mpif_mtime('commtime',0)
924
925! Send variables from getfield process (id_read) to other processes
926!**********************************************************************
927
928! The non-reader processes need to know if cloud water was read.
929    call MPI_Bcast(readclouds,1,MPI_LOGICAL,id_read,MPI_COMM_WORLD,mp_ierr)
930    if (mp_ierr /= 0) goto 600
931
932! Static fields/variables sent only at startup
933    if (first_call) then
934      call MPI_Bcast(oro(:,:),d2_size3,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
935      if (mp_ierr /= 0) goto 600
936      call MPI_Bcast(excessoro(:,:),d2_size3,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
937      if (mp_ierr /= 0) goto 600
938      call MPI_Bcast(lsm(:,:),d2_size3,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
939      if (mp_ierr /= 0) goto 600
940      call MPI_Bcast(xlanduse(:,:,:),d2_size3*numclass,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
941      if (mp_ierr /= 0) goto 600
942      call MPI_Bcast(wftime,d1_size1,MPI_INTEGER,id_read,MPI_COMM_WORLD,mp_ierr)
943      if (mp_ierr /= 0) goto 600
944      call MPI_Bcast(numbwf,1,MPI_INTEGER,id_read,MPI_COMM_WORLD,mp_ierr)
945      if (mp_ierr /= 0) goto 600
946      call MPI_Bcast(nmixz,1,MPI_INTEGER,id_read,MPI_COMM_WORLD,mp_ierr)
947      if (mp_ierr /= 0) goto 600
948      call MPI_Bcast(height,nzmax,MPI_INTEGER,id_read,MPI_COMM_WORLD,mp_ierr)
949      if (mp_ierr /= 0) goto 600
950
951      first_call=.false.
952    endif
953
954    call MPI_Bcast(uu(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
955    if (mp_ierr /= 0) goto 600
956    call MPI_Bcast(vv(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
957    if (mp_ierr /= 0) goto 600
958    call MPI_Bcast(uupol(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
959    if (mp_ierr /= 0) goto 600
960    call MPI_Bcast(vvpol(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
961    if (mp_ierr /= 0) goto 600
962    call MPI_Bcast(ww(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
963    if (mp_ierr /= 0) goto 600
964    call MPI_Bcast(tt(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
965    if (mp_ierr /= 0) goto 600
966    call MPI_Bcast(rho(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
967    if (mp_ierr /= 0) goto 600
968    call MPI_Bcast(drhodz(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
969    if (mp_ierr /= 0) goto 600
970    call MPI_Bcast(tth(:,:,:,li:ui),d3s2,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
971    if (mp_ierr /= 0) goto 600
972    call MPI_Bcast(qvh(:,:,:,li:ui),d3s2,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
973    if (mp_ierr /= 0) goto 600
974    call MPI_Bcast(qv(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
975    if (mp_ierr /= 0) goto 600
976    call MPI_Bcast(pv(:,:,:,li:ui),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
977    if (mp_ierr /= 0) goto 600
978    call MPI_Bcast(clouds(:,:,:,li:ui),d3s1,MPI_INTEGER1,id_read,MPI_COMM_WORLD,mp_ierr)
979    if (mp_ierr /= 0) goto 600
980
981! cloud water/ice:
982    if (readclouds) then
983      call MPI_Bcast(ctwc(:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
984      if (mp_ierr /= 0) goto 600
985    end if
986
987! 2D fields
988    call MPI_Bcast(cloudsh(:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
989    if (mp_ierr /= 0) goto 600
990    call MPI_Bcast(vdep(:,:,:,li:ui),d2s2,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
991    if (mp_ierr /= 0) goto 600
992    call MPI_Bcast(ps(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
993    if (mp_ierr /= 0) goto 600
994    call MPI_Bcast(sd(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
995    if (mp_ierr /= 0) goto 600
996    call MPI_Bcast(tcc(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
997    if (mp_ierr /= 0) goto 600
998    call MPI_Bcast(tt2(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
999    if (mp_ierr /= 0) goto 600
1000    call MPI_Bcast(td2(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1001    if (mp_ierr /= 0) goto 600
1002    call MPI_Bcast(lsprec(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1003    if (mp_ierr /= 0) goto 600
1004    call MPI_Bcast(convprec(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1005    if (mp_ierr /= 0) goto 600
1006    call MPI_Bcast(ustar(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1007    if (mp_ierr /= 0) goto 600
1008    call MPI_Bcast(wstar(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1009    if (mp_ierr /= 0) goto 600
1010    call MPI_Bcast(hmix(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1011    if (mp_ierr /= 0) goto 600
1012    call MPI_Bcast(tropopause(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1013    if (mp_ierr /= 0) goto 600
1014    call MPI_Bcast(oli(:,:,:,li:ui),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1015    if (mp_ierr /= 0) goto 600
1016
1017    if (mp_measure_time) call mpif_mtime('commtime',1)
1018
1019    goto 601
1020
1021600 write(*,*) "mpi_mod> mp_ierr \= 0", mp_ierr
1022    stop
1023
1024601 end subroutine mpif_gf_send_vars
1025
1026
1027  subroutine mpif_gf_send_vars_nest(memstat)
1028!***********************************************************************
1029! DESCRIPTION
1030!   Distribute 'getfield' variables from reader process to all processes.
1031!   For nested fields
1032!
1033!   Called from timemanager
1034!
1035! NOTE
1036!   This subroutine distributes nested windfields from the reader process to
1037!   all other processes. Usually one set of fields is transfered, but at
1038!   first timestep when there are no fields in memory, both are transfered.
1039!   MPI_Bcast is used, so implicitly all processes are synchronized at this
1040!   step
1041!
1042!
1043!***********************************************************************
1044    use com_mod
1045    use par_mod,only: numwfmem
1046
1047    implicit none
1048
1049    integer, intent(in) :: memstat
1050    integer :: i
1051! li,ui    lower,upper indices of fields to be transfered (1-3, ui>=li)
1052    integer :: li,ui
1053
1054! Common array sizes used for communications
1055    integer :: d3_size1 = nxmaxn*nymaxn*nzmax
1056    integer :: d3_size2 = nxmaxn*nymaxn*nuvzmax
1057    integer :: d2_size1 = nxmaxn*nymaxn
1058    integer :: d2_size2 = nxmaxn*nymaxn*maxspec
1059    integer :: d2_size3 = nxmaxn*nymaxn
1060
1061    integer :: d3s1,d3s2,d2s1,d2s2
1062
1063! First time routine is called the unchangeable fields will be transfered   
1064    logical :: first_call=.true.
1065
1066!**********************************************************************
1067
1068! Sizes of arrays transfered
1069    d3s1=d3_size1
1070    d3s2=d3_size2
1071    d2s1=d2_size1
1072    d2s2=d2_size2
1073
1074! Decide which field(s) need to be transfered
1075    if (memstat.ge.32) then ! distribute 2 fields
1076      li=1; ui=2
1077      d3s1=2*d3_size1
1078      d3s2=2*d3_size2
1079      d2s1=2*d2_size1
1080      d2s2=2*d2_size2
1081    else if (memstat.eq.17) then ! distribute 1 field, on 1st index
1082      li=1; ui=1
1083    else if (memstat.eq.18) then ! distribute 1 field, on 2nd index
1084      li=2; ui=2
1085    else if (memstat.eq.19) then ! distribute 1 field, on 3nd index
1086      li=3; ui=3
1087    else
1088      write(*,*) '#### mpi_mod::mpif_gf_send_vars_nest> ERROR: ', &
1089           & 'wrong value of memstat, exiting ####', memstat
1090      stop
1091    end if
1092
1093
1094! Time for MPI communications
1095    if (mp_measure_time) call mpif_mtime('commtime',0)
1096
1097! Distribute variables, send from getfield process (id_read) to other
1098! processes
1099!**********************************************************************
1100
1101! The non-reader processes need to know if cloud water was read.
1102    call MPI_Bcast(readclouds_nest,maxnests,MPI_LOGICAL,id_read,MPI_COMM_WORLD,mp_ierr)
1103    if (mp_ierr /= 0) goto 600
1104
1105! Static fields/variables sent only at startup
1106    if (first_call) then
1107      call MPI_Bcast(oron(:,:,:),d2_size3,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1108      if (mp_ierr /= 0) goto 600
1109      call MPI_Bcast(excessoron(:,:,:),d2_size3,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1110      if (mp_ierr /= 0) goto 600
1111      call MPI_Bcast(lsmn(:,:,:),d2_size3,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1112      if (mp_ierr /= 0) goto 600
1113      call MPI_Bcast(xlandusen(:,:,:,:),d2_size3*numclass,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1114      if (mp_ierr /= 0) goto 600
1115      first_call=.false.
1116    end if
1117
1118! MPI prefers contiguous arrays for sending (else a buffer is created),
1119! hence the loop over nests
1120!**********************************************************************
1121    do i=1, numbnests
1122! 3D fields
1123      call MPI_Bcast(uun(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1124      if (mp_ierr /= 0) goto 600
1125      call MPI_Bcast(vvn(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1126      if (mp_ierr /= 0) goto 600
1127      call MPI_Bcast(wwn(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1128      if (mp_ierr /= 0) goto 600
1129      call MPI_Bcast(ttn(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1130      if (mp_ierr /= 0) goto 600
1131      call MPI_Bcast(rhon(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1132      if (mp_ierr /= 0) goto 600
1133      call MPI_Bcast(drhodzn(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1134      if (mp_ierr /= 0) goto 600
1135      call MPI_Bcast(tthn(:,:,:,li:ui,i),d3s2,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1136      if (mp_ierr /= 0) goto 600
1137      call MPI_Bcast(qvhn(:,:,:,li:ui,i),d3s2,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1138      if (mp_ierr /= 0) goto 600
1139      call MPI_Bcast(qvn(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1140      if (mp_ierr /= 0) goto 600
1141      call MPI_Bcast(pvn(:,:,:,li:ui,i),d3s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1142      if (mp_ierr /= 0) goto 600
1143      call MPI_Bcast(cloudsn(:,:,:,li:ui,i),d3s1,MPI_INTEGER1,id_read,MPI_COMM_WORLD,mp_ierr)
1144      if (mp_ierr /= 0) goto 600
1145
1146! cloud water/ice:
1147    if (readclouds_nest(i)) then
1148      ! call MPI_Bcast(icloud_stats(:,:,:,li:ui),d2s1*5,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1149      ! if (mp_ierr /= 0) goto 600
1150      call MPI_Bcast(ctwcn(:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1151      if (mp_ierr /= 0) goto 600
1152    end if
1153
1154! 2D fields
1155      call MPI_Bcast(cloudshn(:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1156      if (mp_ierr /= 0) goto 600
1157      call MPI_Bcast(vdepn(:,:,:,li:ui,i),d2s2,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1158      if (mp_ierr /= 0) goto 600
1159      call MPI_Bcast(psn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1160      if (mp_ierr /= 0) goto 600
1161      call MPI_Bcast(sdn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1162      if (mp_ierr /= 0) goto 600
1163      call MPI_Bcast(tccn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1164      if (mp_ierr /= 0) goto 600
1165      call MPI_Bcast(tt2n(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1166      if (mp_ierr /= 0) goto 600
1167      call MPI_Bcast(td2n(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1168      if (mp_ierr /= 0) goto 600
1169      call MPI_Bcast(lsprecn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1170      if (mp_ierr /= 0) goto 600
1171      call MPI_Bcast(convprecn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1172      if (mp_ierr /= 0) goto 600
1173      call MPI_Bcast(ustarn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1174      if (mp_ierr /= 0) goto 600
1175      call MPI_Bcast(wstarn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1176      if (mp_ierr /= 0) goto 600
1177      call MPI_Bcast(olin(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1178      if (mp_ierr /= 0) goto 600
1179      call MPI_Bcast(hmixn(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1180      if (mp_ierr /= 0) goto 600
1181      call MPI_Bcast(tropopausen(:,:,:,li:ui,i),d2s1,mp_sp,id_read,MPI_COMM_WORLD,mp_ierr)
1182      if (mp_ierr /= 0) goto 600
1183    end do
1184
1185    if (mp_measure_time) call mpif_mtime('commtime',1)
1186
1187    goto 601
1188
1189600 write(*,*) "mpi_mod> mp_ierr \= 0",mp_ierr
1190    stop
1191
1192601 end subroutine mpif_gf_send_vars_nest
1193
1194
1195  subroutine mpif_gf_send_vars_async(memstat)
1196!*******************************************************************************
1197! DESCRIPTION
1198!   Distribute 'getfield' variables from reader process to all processes.
1199!   Called from timemanager by reader process only
1200!
1201! NOTE
1202!   This version uses asynchronious sends. The newest fields are sent in the
1203!   background, so calculations can continue while
1204!   MPI communications are performed.
1205!
1206!   The way the MPI tags/requests are sequenced means that this routine must
1207!   carefully match routine 'mpif_gf_recv_vars_async'
1208!
1209! VARIABLES
1210!   memstat -- input, for resolving pointer to windfield index being read
1211!   mind    -- index where to place new fields
1212!
1213!
1214!*******************************************************************************
1215    use com_mod
1216
1217    implicit none
1218
1219    integer, intent(in) :: memstat
1220    integer :: mind
1221    integer :: dest,i
1222
1223! Common array sizes used for communications
1224    integer :: d3s1 = nxmax*nymax*nzmax
1225    integer :: d3s2 = nxmax*nymax*nuvzmax
1226    integer :: d2s1 = nxmax*nymax
1227    integer :: d2s2 = nxmax*nymax*maxspec
1228!    integer :: d1s1 = maxwf
1229
1230!*******************************************************************************
1231
1232! At the time the send is posted, the reader process is one step further
1233! in the permutation of memstat compared with the receiving processes
1234
1235    if (memstat.ge.32) then
1236! last read was synchronous, to indices 1 and 2, 3 is free
1237      write(*,*) "#### mpi_mod::mpif_gf_send_vars_async> ERROR: &
1238           & memstat>=32 should never happen here."
1239      stop
1240    else if (memstat.eq.17) then
1241! old fields on 1,2, send 3
1242      mind=3
1243    else if (memstat.eq.18) then
1244! old fields on 2,3, send 1
1245      mind=1
1246    else if (memstat.eq.19) then
1247! old fields on 3,1, send 2
1248      mind=2
1249    else
1250      write(*,*) "#### mpi_mod::mpif_gf_send_vars_async> ERROR: &
1251           & invalid memstat"
1252    end if
1253
1254    if (mp_dev_mode) then
1255      if (mp_pid.ne.id_read) then
1256        write(*,*) 'MPI_DEV_MODE: error calling mpif_gf_send_vars_async'
1257      end if
1258    end if
1259
1260    if (mp_dev_mode) write(*,*) '## in mpif_gf_send_vars_async, memstat:', memstat
1261
1262! Time for MPI communications
1263    if (mp_measure_time) call mpif_mtime('commtime',0)
1264
1265! Loop over receiving processes, initiate data sending
1266!*****************************************************
1267
1268    do dest=0,mp_np-1 ! mp_np-2 will also work if last proc reserved for reading
1269! TODO: use mp_partgroup_np here
1270      if (dest.eq.id_read) cycle
1271      i=dest*nvar_async
1272      call MPI_Isend(uu(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1273      if (mp_ierr /= 0) goto 600
1274      i=i+1
1275      call MPI_Isend(vv(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1276      if (mp_ierr /= 0) goto 600
1277      i=i+1
1278      call MPI_Isend(uupol(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1279      if (mp_ierr /= 0) goto 600
1280      i=i+1
1281      call MPI_Isend(vvpol(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1282      if (mp_ierr /= 0) goto 600
1283      i=i+1
1284      call MPI_Isend(ww(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1285      if (mp_ierr /= 0) goto 600
1286      i=i+1
1287      call MPI_Isend(tt(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1288      if (mp_ierr /= 0) goto 600
1289      i=i+1
1290      call MPI_Isend(rho(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1291      if (mp_ierr /= 0) goto 600
1292      i=i+1
1293      call MPI_Isend(drhodz(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1294      if (mp_ierr /= 0) goto 600
1295      i=i+1
1296      call MPI_Isend(tth(:,:,:,mind),d3s2,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1297      if (mp_ierr /= 0) goto 600
1298      i=i+1
1299      call MPI_Isend(qvh(:,:,:,mind),d3s2,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1300      if (mp_ierr /= 0) goto 600
1301      i=i+1
1302      call MPI_Isend(qv(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1303      if (mp_ierr /= 0) goto 600
1304      i=i+1
1305      call MPI_Isend(pv(:,:,:,mind),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1306      if (mp_ierr /= 0) goto 600
1307      i=i+1
1308      call MPI_Isend(clouds(:,:,:,mind),d3s1,MPI_INTEGER1,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1309      i=i+1
1310      if (mp_ierr /= 0) goto 600
1311      call MPI_Isend(cloudsh(:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1312      if (mp_ierr /= 0) goto 600
1313      i=i+1
1314      call MPI_Isend(vdep(:,:,:,mind),d2s2,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1315      if (mp_ierr /= 0) goto 600
1316      i=i+1
1317! 15
1318      call MPI_Isend(ps(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1319      if (mp_ierr /= 0) goto 600
1320      i=i+1
1321      call MPI_Isend(sd(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1322      if (mp_ierr /= 0) goto 600
1323      i=i+1
1324      call MPI_Isend(tcc(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1325      if (mp_ierr /= 0) goto 600
1326      i=i+1
1327      call MPI_Isend(tt2(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1328      if (mp_ierr /= 0) goto 600
1329      i=i+1
1330      call MPI_Isend(td2(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1331      if (mp_ierr /= 0) goto 600
1332      i=i+1
1333      call MPI_Isend(lsprec(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1334      if (mp_ierr /= 0) goto 600
1335      i=i+1
1336      call MPI_Isend(convprec(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1337      if (mp_ierr /= 0) goto 600
1338      i=i+1
1339      call MPI_Isend(ustar(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1340      if (mp_ierr /= 0) goto 600
1341      i=i+1
1342      call MPI_Isend(wstar(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1343      if (mp_ierr /= 0) goto 600
1344      i=i+1
1345      call MPI_Isend(hmix(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1346      if (mp_ierr /= 0) goto 600
1347      i=i+1
1348! 25
1349      call MPI_Isend(tropopause(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1350      if (mp_ierr /= 0) goto 600
1351      i=i+1
1352      call MPI_Isend(oli(:,:,:,mind),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1353      if (mp_ierr /= 0) goto 600
1354
1355! Send cloud water if it exists. Increment counter always (as on receiving end)
1356      if (readclouds) then
1357        i=i+1
1358        call MPI_Isend(ctwc(:,:,mind),d2s1,mp_sp,dest,tm1,&
1359             &MPI_COMM_WORLD,reqs(i),mp_ierr)
1360        if (mp_ierr /= 0) goto 600
1361      end if
1362    end do
1363
1364    if (mp_measure_time) call mpif_mtime('commtime',1)
1365
1366    goto 601
1367
1368600 write(*,*) "#### mpi_mod::mpif_gf_send_vars_async> mp_ierr \= 0", mp_ierr
1369    stop
1370
1371601 end subroutine mpif_gf_send_vars_async
1372
1373
1374  subroutine mpif_gf_recv_vars_async(memstat)
1375!*******************************************************************************
1376! DESCRIPTION
1377!   Receive 'getfield' variables from reader process.
1378!   Called from timemanager by all processes except reader
1379!
1380! NOTE
1381!   This version uses asynchronious communications.
1382!
1383! VARIABLES
1384!   memstat -- input, used to resolve windfield index being received
1385!
1386!
1387!*******************************************************************************
1388    use com_mod
1389
1390    implicit none
1391
1392    integer, intent(in) :: memstat
1393    integer :: mind,j
1394
1395! Common array sizes used for communications
1396    integer :: d3s1 = nxmax*nymax*nzmax
1397    integer :: d3s2 = nxmax*nymax*nuvzmax
1398    integer :: d2s1 = nxmax*nymax
1399    integer :: d2s2 = nxmax*nymax*maxspec
1400    !integer :: d1_size1 = maxwf
1401
1402!    integer :: d3s1,d3s2,d2s1,d2s2
1403!*******************************************************************************
1404
1405! At the time this immediate receive is posted, memstat is the state of
1406! windfield indices at the previous/current time. From this, the future
1407! state is deduced.
1408    if (memstat.eq.32) then
1409! last read was synchronous to indices 1 and 2, 3 is free
1410      mind=3
1411    else if (memstat.eq.17) then
1412! last read was asynchronous to index 3, 1 is free
1413      mind=1
1414    else if (memstat.eq.18) then
1415! last read was asynchronous to index 1, 2 is free
1416      mind=2
1417    else if (memstat.eq.19) then
1418! last read was asynchronous to index 2, 3 is free
1419      mind=3
1420    else
1421! illegal state
1422      write(*,*) 'mpi_mod> FLEXPART ERROR: Illegal memstat value. Exiting.'
1423      stop
1424    end if
1425
1426    if (mp_dev_mode) then
1427      if (mp_pid.eq.id_read) then
1428        write(*,*) 'MPI_DEV_MODE: error calling mpif_gf_recv_vars_async'
1429      end if
1430    end if
1431
1432! Time for MPI communications
1433    if (mp_measure_time) call mpif_mtime('commtime',0)
1434
1435    if (mp_dev_mode) write(*,*) '## in mpif_gf_send_vars_async, memstat:', memstat
1436
1437! Initiate receiving of data
1438!***************************
1439
1440! Get MPI tags/requests for communications
1441    j=mp_pid*nvar_async
1442    call MPI_Irecv(uu(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1443         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1444    if (mp_ierr /= 0) goto 600
1445    j=j+1
1446    call MPI_Irecv(vv(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1447         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1448    if (mp_ierr /= 0) goto 600
1449    j=j+1
1450    call MPI_Irecv(uupol(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1451         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1452    if (mp_ierr /= 0) goto 600
1453    j=j+1
1454    call MPI_Irecv(vvpol(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1455         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1456    if (mp_ierr /= 0) goto 600
1457    j=j+1
1458    call MPI_Irecv(ww(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1459         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1460    if (mp_ierr /= 0) goto 600
1461    j=j+1
1462    call MPI_Irecv(tt(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1463         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1464    if (mp_ierr /= 0) goto 600
1465    j=j+1
1466    call MPI_Irecv(rho(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1467         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1468    if (mp_ierr /= 0) goto 600
1469    j=j+1
1470    call MPI_Irecv(drhodz(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1471         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1472    if (mp_ierr /= 0) goto 600
1473    j=j+1
1474    call MPI_Irecv(tth(:,:,:,mind),d3s2,mp_sp,id_read,MPI_ANY_TAG,&
1475         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1476    if (mp_ierr /= 0) goto 600
1477    j=j+1
1478    call MPI_Irecv(qvh(:,:,:,mind),d3s2,mp_sp,id_read,MPI_ANY_TAG,&
1479         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1480    if (mp_ierr /= 0) goto 600
1481    j=j+1
1482    call MPI_Irecv(qv(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1483         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1484    if (mp_ierr /= 0) goto 600
1485    j=j+1
1486    call MPI_Irecv(pv(:,:,:,mind),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1487         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1488    if (mp_ierr /= 0) goto 600
1489    j=j+1
1490    call MPI_Irecv(clouds(:,:,:,mind),d3s1,MPI_INTEGER1,id_read,MPI_ANY_TAG,&
1491         &MPI_COMM_WORLD,reqs(j),mp_ierr)   
1492    if (mp_ierr /= 0) goto 600
1493    j=j+1
1494    call MPI_Irecv(cloudsh(:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1495         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1496    if (mp_ierr /= 0) goto 600
1497    j=j+1
1498    call MPI_Irecv(vdep(:,:,:,mind),d2s2,mp_sp,id_read,MPI_ANY_TAG,&
1499         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1500    if (mp_ierr /= 0) goto 600
1501    j=j+1
1502    call MPI_Irecv(ps(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1503         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1504    if (mp_ierr /= 0) goto 600
1505    j=j+1
1506    call MPI_Irecv(sd(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1507         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1508    if (mp_ierr /= 0) goto 600
1509    j=j+1
1510    call MPI_Irecv(tcc(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1511         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1512    if (mp_ierr /= 0) goto 600
1513    j=j+1
1514    call MPI_Irecv(tt2(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1515         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1516    if (mp_ierr /= 0) goto 600
1517    j=j+1
1518    call MPI_Irecv(td2(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1519         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1520    if (mp_ierr /= 0) goto 600
1521    j=j+1
1522    call MPI_Irecv(lsprec(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1523         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1524    if (mp_ierr /= 0) goto 600
1525    j=j+1
1526    call MPI_Irecv(convprec(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1527         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1528    if (mp_ierr /= 0) goto 600
1529    call MPI_Irecv(ustar(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1530         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1531    if (mp_ierr /= 0) goto 600
1532    j=j+1
1533    call MPI_Irecv(wstar(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1534         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1535    if (mp_ierr /= 0) goto 600
1536    j=j+1
1537    call MPI_Irecv(hmix(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1538         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1539    if (mp_ierr /= 0) goto 600
1540    j=j+1
1541    call MPI_Irecv(tropopause(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1542         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1543    if (mp_ierr /= 0) goto 600
1544    j=j+1
1545    call MPI_Irecv(oli(:,:,:,mind),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1546         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1547    if (mp_ierr /= 0) goto 600
1548
1549! Post request for clwc. These data are possibly not sent, request must then be cancelled
1550! For now assume that data at all steps either have or do not have water
1551    if (readclouds) then
1552      j=j+1
1553      call MPI_Irecv(ctwc(:,:,mind),d2s1*5,mp_sp,id_read,MPI_ANY_TAG,&
1554           &MPI_COMM_WORLD,reqs(j),mp_ierr)
1555      if (mp_ierr /= 0) goto 600
1556    end if
1557
1558
1559    if (mp_measure_time) call mpif_mtime('commtime',1)
1560
1561    goto 601
1562
1563600 write(*,*) "#### mpi_mod::mpif_gf_recv_vars_async> MPI ERROR ", mp_ierr
1564    stop
1565
1566601 end subroutine mpif_gf_recv_vars_async
1567
1568
1569  subroutine mpif_gf_send_vars_nest_async(memstat)
1570!*******************************************************************************
1571! DESCRIPTION
1572!   Distribute 'getfield' variables from reader process to all processes.
1573!   Called from timemanager by reader process only.
1574!   Version for nested wind fields
1575!
1576! NOTE
1577!   This version uses asynchronious sends. The newest fields are sent in the
1578!   background, so calculations can continue while
1579!   MPI communications are performed.
1580!
1581!   The way the MPI tags/requests are sequenced means that this routine must
1582!   carefully match routine 'mpif_gf_recv_vars_async'
1583!
1584! VARIABLES
1585!   memstat -- input, for resolving pointer to windfield index being read
1586!   mind    -- index where to place new fields
1587!
1588! TODO
1589!   Some unused arrays are currently sent (uupoln,..)
1590!*******************************************************************************
1591    use com_mod
1592
1593    implicit none
1594
1595    integer, intent(in) :: memstat
1596    integer :: mind
1597    integer :: dest,i,k
1598
1599! Common array sizes used for communications
1600    integer :: d3s1 = nxmaxn*nymaxn*nzmax
1601    integer :: d3s2 = nxmaxn*nymaxn*nuvzmax
1602    integer :: d2s1 = nxmaxn*nymaxn
1603    integer :: d2s2 = nxmaxn*nymaxn*maxspec
1604
1605!*******************************************************************************
1606
1607! At the time the send is posted, the reader process is one step further
1608! in the permutation of memstat compared with the receiving processes
1609
1610    if (memstat.ge.32) then
1611! last read was synchronous, to indices 1 and 2, 3 is free
1612      write(*,*) "#### mpi_mod::mpif_gf_send_vars_nest_async> ERROR: &
1613           & memstat>=32 should never happen here."
1614      stop
1615    else if (memstat.eq.17) then
1616! old fields on 1,2, send 3
1617      mind=3
1618    else if (memstat.eq.18) then
1619! old fields on 2,3, send 1
1620      mind=1
1621    else if (memstat.eq.19) then
1622! old fields on 3,1, send 2
1623      mind=2
1624    else
1625      write(*,*) "#### mpi_mod::mpif_gf_send_vars_nest_async> ERROR: &
1626           & invalid memstat"
1627    end if
1628
1629    if (mp_dev_mode) then
1630      if (mp_pid.ne.id_read) then
1631        write(*,*) 'MPI_DEV_MODE: error calling mpif_gf_send_vars_nest_async'
1632      end if
1633    end if
1634
1635    if (mp_dev_mode) write(*,*) '## in mpif_gf_send_vars_nest_async, memstat:', memstat
1636
1637! Time for MPI communications
1638    if (mp_measure_time) call mpif_mtime('commtime',0)
1639
1640! Loop over receiving processes, initiate data sending
1641!*****************************************************
1642
1643    do dest=0,mp_np-1 ! mp_np-2 will also work if last proc reserved for reading
1644! TODO: use mp_partgroup_np here
1645      if (dest.eq.id_read) cycle
1646      do k=1, numbnests
1647      i=dest*nvar_async
1648      call MPI_Isend(uun(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1649      if (mp_ierr /= 0) goto 600
1650      i=i+1
1651      call MPI_Isend(vvn(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1652      if (mp_ierr /= 0) goto 600
1653      i=i+1
1654      call MPI_Isend(wwn(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1655      if (mp_ierr /= 0) goto 600
1656      i=i+1
1657      call MPI_Isend(ttn(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1658      if (mp_ierr /= 0) goto 600
1659      i=i+1
1660      call MPI_Isend(rhon(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1661      if (mp_ierr /= 0) goto 600
1662      i=i+1
1663      call MPI_Isend(drhodzn(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1664      if (mp_ierr /= 0) goto 600
1665      i=i+1
1666      call MPI_Isend(tthn(:,:,:,mind,k),d3s2,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1667      if (mp_ierr /= 0) goto 600
1668      i=i+1
1669      call MPI_Isend(qvhn(:,:,:,mind,k),d3s2,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1670      if (mp_ierr /= 0) goto 600
1671      i=i+1
1672      call MPI_Isend(qvn(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1673      if (mp_ierr /= 0) goto 600
1674      i=i+1
1675      call MPI_Isend(pvn(:,:,:,mind,k),d3s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1676      if (mp_ierr /= 0) goto 600
1677      i=i+1
1678      call MPI_Isend(cloudsn(:,:,:,mind,k),d3s1,MPI_INTEGER1,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1679      i=i+1
1680      if (mp_ierr /= 0) goto 600
1681      call MPI_Isend(cloudshn(:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1682      if (mp_ierr /= 0) goto 600
1683      i=i+1
1684      call MPI_Isend(vdepn(:,:,:,mind,k),d2s2,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1685      if (mp_ierr /= 0) goto 600
1686      i=i+1
1687      call MPI_Isend(psn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1688      if (mp_ierr /= 0) goto 600
1689      i=i+1
1690      call MPI_Isend(sdn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1691      if (mp_ierr /= 0) goto 600
1692      i=i+1
1693! 15
1694      call MPI_Isend(tccn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1695      if (mp_ierr /= 0) goto 600
1696      i=i+1
1697      call MPI_Isend(tt2n(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1698      if (mp_ierr /= 0) goto 600
1699      i=i+1
1700      call MPI_Isend(td2n(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1701      if (mp_ierr /= 0) goto 600
1702      i=i+1
1703      call MPI_Isend(lsprecn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1704      if (mp_ierr /= 0) goto 600
1705      i=i+1
1706      call MPI_Isend(convprecn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1707      if (mp_ierr /= 0) goto 600
1708      i=i+1
1709      call MPI_Isend(ustarn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1710      if (mp_ierr /= 0) goto 600
1711      i=i+1
1712      call MPI_Isend(wstarn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1713      if (mp_ierr /= 0) goto 600
1714      i=i+1
1715      call MPI_Isend(hmixn(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1716      if (mp_ierr /= 0) goto 600
1717      i=i+1
1718      call MPI_Isend(tropopausen(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1719      if (mp_ierr /= 0) goto 600
1720      i=i+1
1721      call MPI_Isend(olin(:,:,:,mind,k),d2s1,mp_sp,dest,tm1,MPI_COMM_WORLD,reqs(i),mp_ierr)
1722      if (mp_ierr /= 0) goto 600
1723! 25
1724
1725! Send cloud water if it exists. Increment counter always (as on receiving end)
1726      if (readclouds) then
1727        i=i+1
1728        call MPI_Isend(ctwcn(:,:,mind,k),d2s1,mp_sp,dest,tm1,&
1729             &MPI_COMM_WORLD,reqs(i),mp_ierr)
1730        if (mp_ierr /= 0) goto 600
1731      end if
1732    end do
1733  end do
1734
1735    if (mp_measure_time) call mpif_mtime('commtime',1)
1736
1737    goto 601
1738
1739600 write(*,*) "#### mpi_mod::mpif_gf_send_vars_nest_async> mp_ierr \= 0", mp_ierr
1740    stop
1741
1742601 end subroutine mpif_gf_send_vars_nest_async
1743
1744
1745  subroutine mpif_gf_recv_vars_nest_async(memstat)
1746!*******************************************************************************
1747! DESCRIPTION
1748!   Receive 'getfield' variables from reader process.
1749!   Called from timemanager by all processes except reader
1750!   Version for nested wind fields
1751!
1752! NOTE
1753!   This version uses asynchronious communications.
1754!
1755! VARIABLES
1756!   memstat -- input, used to resolve windfield index being received
1757!
1758!
1759!*******************************************************************************
1760    use com_mod
1761
1762    implicit none
1763
1764    integer, intent(in) :: memstat
1765    integer :: mind,j,k
1766
1767! Common array sizes used for communications
1768    integer :: d3s1 = nxmaxn*nymaxn*nzmax
1769    integer :: d3s2 = nxmaxn*nymaxn*nuvzmax
1770    integer :: d2s1 = nxmaxn*nymaxn
1771    integer :: d2s2 = nxmaxn*nymaxn*maxspec
1772
1773!*******************************************************************************
1774
1775! At the time this immediate receive is posted, memstat is the state of
1776! windfield indices at the previous/current time. From this, the future
1777! state is deduced.
1778    if (memstat.eq.32) then
1779! last read was synchronous to indices 1 and 2, 3 is free
1780      mind=3
1781    else if (memstat.eq.17) then
1782! last read was asynchronous to index 3, 1 is free
1783      mind=1
1784    else if (memstat.eq.18) then
1785! last read was asynchronous to index 1, 2 is free
1786      mind=2
1787    else if (memstat.eq.19) then
1788! last read was asynchronous to index 2, 3 is free
1789      mind=3
1790    else
1791! illegal state
1792      write(*,*) 'mpi_mod> FLEXPART ERROR: Illegal memstat value. Exiting.'
1793      stop
1794    end if
1795
1796    if (mp_dev_mode) then
1797      if (mp_pid.eq.id_read) then
1798        write(*,*) 'MPI_DEV_MODE: error calling mpif_gf_recv_vars_async'
1799      end if
1800    end if
1801
1802! Time for MPI communications
1803    if (mp_measure_time) call mpif_mtime('commtime',0)
1804
1805    if (mp_dev_mode) write(*,*) '## in mpif_gf_send_vars_async, memstat:', memstat
1806
1807! Initiate receiving of data
1808!***************************
1809
1810    do k=1, numbnests
1811! Get MPI tags/requests for communications
1812    j=mp_pid*nvar_async
1813    call MPI_Irecv(uun(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1814         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1815    if (mp_ierr /= 0) goto 600
1816    j=j+1
1817    call MPI_Irecv(vvn(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1818         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1819    if (mp_ierr /= 0) goto 600
1820    j=j+1
1821    call MPI_Irecv(wwn(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1822         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1823    if (mp_ierr /= 0) goto 600
1824    j=j+1
1825    call MPI_Irecv(ttn(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1826         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1827    if (mp_ierr /= 0) goto 600
1828    j=j+1
1829    call MPI_Irecv(rhon(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1830         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1831    if (mp_ierr /= 0) goto 600
1832    j=j+1
1833    call MPI_Irecv(drhodzn(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1834         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1835    if (mp_ierr /= 0) goto 600
1836    j=j+1
1837    call MPI_Irecv(tthn(:,:,:,mind,k),d3s2,mp_sp,id_read,MPI_ANY_TAG,&
1838         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1839    if (mp_ierr /= 0) goto 600
1840    j=j+1
1841    call MPI_Irecv(qvhn(:,:,:,mind,k),d3s2,mp_sp,id_read,MPI_ANY_TAG,&
1842         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1843    if (mp_ierr /= 0) goto 600
1844    j=j+1
1845    call MPI_Irecv(qvn(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1846         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1847    if (mp_ierr /= 0) goto 600
1848    j=j+1
1849    call MPI_Irecv(pvn(:,:,:,mind,k),d3s1,mp_sp,id_read,MPI_ANY_TAG,&
1850         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1851    if (mp_ierr /= 0) goto 600
1852    j=j+1
1853    call MPI_Irecv(cloudsn(:,:,:,mind,k),d3s1,MPI_INTEGER1,id_read,MPI_ANY_TAG,&
1854         &MPI_COMM_WORLD,reqs(j),mp_ierr)   
1855    if (mp_ierr /= 0) goto 600
1856    j=j+1
1857    call MPI_Irecv(cloudshn(:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1858         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1859    if (mp_ierr /= 0) goto 600
1860    j=j+1
1861    call MPI_Irecv(vdepn(:,:,:,mind,k),d2s2,mp_sp,id_read,MPI_ANY_TAG,&
1862         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1863    if (mp_ierr /= 0) goto 600
1864    j=j+1
1865    call MPI_Irecv(psn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1866         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1867    if (mp_ierr /= 0) goto 600
1868    j=j+1
1869    call MPI_Irecv(sdn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1870         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1871    if (mp_ierr /= 0) goto 600
1872    j=j+1
1873    call MPI_Irecv(tccn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1874         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1875    if (mp_ierr /= 0) goto 600
1876    j=j+1
1877    call MPI_Irecv(tt2n(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1878         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1879    if (mp_ierr /= 0) goto 600
1880    j=j+1
1881    call MPI_Irecv(td2n(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1882         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1883    if (mp_ierr /= 0) goto 600
1884    j=j+1
1885    call MPI_Irecv(lsprecn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1886         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1887    if (mp_ierr /= 0) goto 600
1888    j=j+1
1889    call MPI_Irecv(convprecn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1890         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1891    if (mp_ierr /= 0) goto 600
1892    call MPI_Irecv(ustarn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1893         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1894    if (mp_ierr /= 0) goto 600
1895    j=j+1
1896    call MPI_Irecv(wstarn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1897         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1898    if (mp_ierr /= 0) goto 600
1899    j=j+1
1900    call MPI_Irecv(hmixn(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1901         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1902    if (mp_ierr /= 0) goto 600
1903    j=j+1
1904    call MPI_Irecv(tropopausen(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1905         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1906    if (mp_ierr /= 0) goto 600
1907    j=j+1
1908    call MPI_Irecv(olin(:,:,:,mind,k),d2s1,mp_sp,id_read,MPI_ANY_TAG,&
1909         &MPI_COMM_WORLD,reqs(j),mp_ierr)
1910    if (mp_ierr /= 0) goto 600
1911
1912! Post request for clwc. These data are possibly not sent, request must then be cancelled
1913! For now assume that data at all steps either have or do not have water
1914    if (readclouds) then
1915      j=j+1
1916      call MPI_Irecv(ctwcn(:,:,mind,k),d2s1*5,mp_sp,id_read,MPI_ANY_TAG,&
1917           &MPI_COMM_WORLD,reqs(j),mp_ierr)
1918      if (mp_ierr /= 0) goto 600
1919    end if
1920  end do
1921
1922    if (mp_measure_time) call mpif_mtime('commtime',1)
1923
1924    goto 601
1925
1926600 write(*,*) "#### mpi_mod::mpif_gf_recv_vars_nest_async> MPI ERROR ", mp_ierr
1927    stop
1928
1929601 end subroutine mpif_gf_recv_vars_nest_async
1930
1931
1932  subroutine mpif_gf_request
1933!*******************************************************************************
1934! DESCRIPTION
1935!   Check for completion of MPI_Isend/MPI_Irecv data transfer.
1936!   
1937!   
1938! NOTE
1939!   implicit synchronisation between all processes takes place here
1940!
1941! TODO
1942!   take into account nested wind fields
1943!
1944!*******************************************************************************
1945!    use com_mod,only: readclouds
1946
1947    implicit none
1948
1949
1950    integer :: n_req !,i
1951
1952!***********************************************************************
1953
1954    n_req=nvar_async*mp_np
1955
1956    if (mp_measure_time) call mpif_mtime('commtime',0)
1957
1958!    call MPI_Wait(rm1,MPI_STATUS_IGNORE,mp_ierr)
1959
1960! TODO: cancel recv request if readclouds=.false.
1961!    if (readclouds) then
1962    call MPI_Waitall(n_req,reqs,MPI_STATUSES_IGNORE,mp_ierr)
1963!    endif
1964! else
1965!   do i = 0, nvar_async*mp_np-1
1966!     if (mod(i,27).eq.0 .or. mod(i,28).eq.0) then
1967!       call MPI_Cancel(reqs(i),mp_ierr)
1968!       cycle
1969!     end if
1970!     call MPI_Wait(reqs(i),MPI_STATUS_IGNORE,mp_ierr)
1971!   end do
1972! end if
1973
1974    if (mp_ierr /= 0) goto 600
1975
1976    if (mp_measure_time) call mpif_mtime('commtime',1)
1977
1978    goto 601
1979
1980600 write(*,*) "#### mpi_mod::mpif_gf_request> MPI ERROR ", mp_ierr
1981    stop
1982
1983601 end subroutine mpif_gf_request
1984
1985
1986  subroutine mpif_tm_reduce_grid
1987!***********************************************************************
1988! Collect grid variable to PID 0, adding from all processes.
1989!
1990! NOTE
1991!   - Older MPI implementations may lack the MPI_IN_PLACE operation.
1992!     If so use 1) below and change unc_mod
1993!
1994!***********************************************************************
1995    use com_mod
1996    use unc_mod
1997    use par_mod
1998
1999    implicit none
2000
2001    integer :: grid_size2d,grid_size3d
2002    integer, parameter :: rcpt_size=maxreceptor*maxspec
2003
2004!**********************************************************************
2005    grid_size2d=numxgrid*numygrid*maxspec* &
2006         & maxpointspec_act*nclassunc*maxageclass
2007    grid_size3d=numxgrid*numygrid*numzgrid*maxspec* &
2008         & maxpointspec_act*nclassunc*maxageclass
2009
2010
2011! Time for MPI communications
2012    if (mp_measure_time) call mpif_mtime('commtime',0)
2013
2014! 1) Using a separate grid (gridunc0) for received values
2015! call MPI_Reduce(gridunc, gridunc0, grid_size3d, mp_sp, MPI_SUM, id_root, &
2016!      & mp_comm_used, mp_ierr)
2017! if (mp_ierr /= 0) goto 600
2018
2019! 2) Using in-place reduction
2020    if (lroot) then
2021      call MPI_Reduce(MPI_IN_PLACE, gridunc, grid_size3d, mp_sp, MPI_SUM, id_root, &
2022           & mp_comm_used, mp_ierr)
2023      if (mp_ierr /= 0) goto 600
2024    else
2025      call MPI_Reduce(gridunc, 0, grid_size3d, mp_sp, MPI_SUM, id_root, &
2026           & mp_comm_used, mp_ierr)
2027    end if
2028
2029    if ((WETDEP).and.(ldirect.gt.0)) then
2030      call MPI_Reduce(wetgridunc, wetgridunc0, grid_size2d, mp_cp, MPI_SUM, id_root, &
2031           & mp_comm_used, mp_ierr)
2032      if (mp_ierr /= 0) goto 600
2033    end if
2034
2035    if ((DRYDEP).and.(ldirect.gt.0)) then
2036      call MPI_Reduce(drygridunc, drygridunc0, grid_size2d, mp_cp, MPI_SUM, id_root, &
2037           & mp_comm_used, mp_ierr)
2038      if (mp_ierr /= 0) goto 600
2039    end if
2040
2041! Receptor concentrations   
2042    if (lroot) then
2043      call MPI_Reduce(MPI_IN_PLACE,creceptor,rcpt_size,mp_sp,MPI_SUM,id_root, &
2044           & mp_comm_used,mp_ierr)
2045      if (mp_ierr /= 0) goto 600
2046    else
2047      call MPI_Reduce(creceptor,0,rcpt_size,mp_sp,MPI_SUM,id_root, &
2048           & mp_comm_used,mp_ierr)
2049    end if
2050
2051    if (mp_measure_time) call mpif_mtime('commtime',1)
2052
2053    goto 601
2054
2055600 write(*,*) "mpi_mod> mp_ierr \= 0", mp_ierr
2056    stop
2057
2058601 end subroutine mpif_tm_reduce_grid
2059
2060
2061  subroutine mpif_tm_reduce_grid_nest
2062!***********************************************************************
2063! Collect nested grids to PID 0, adding from all processes.
2064!
2065! NOTE
2066!   - Compiling with 'fcheck=all' (gfortran) will cause a run-time error
2067!     as wetgriduncn0 is only allocated for root process. Possibly valid
2068!     MPI code but invalid gfortran.
2069!
2070!***********************************************************************
2071    use com_mod
2072    use unc_mod
2073    use par_mod
2074
2075    implicit none
2076
2077    integer :: grid_size2d,grid_size3d
2078
2079!**********************************************************************
2080
2081    grid_size3d=numxgridn*numygridn*numzgrid*maxspec* &
2082         & maxpointspec_act*nclassunc*maxageclass
2083    grid_size2d=numxgridn*numygridn*maxspec* &
2084         & maxpointspec_act*nclassunc*maxageclass
2085
2086
2087! Time for MPI communications
2088    if (mp_measure_time) call mpif_mtime('commtime',0)
2089
2090! Using a separate grid (gridunc0) for received values, for debugging
2091! call MPI_Reduce(griduncn, griduncn0, grid_size3d, mp_sp, MPI_SUM, id_root, &
2092!      & mp_comm_used, mp_ierr)
2093! if (mp_ierr /= 0) goto 600
2094
2095! Using in-place reduction
2096    if (lroot) then
2097      call MPI_Reduce(MPI_IN_PLACE, griduncn, grid_size3d, mp_sp, MPI_SUM, id_root, &
2098           & mp_comm_used, mp_ierr)
2099      if (mp_ierr /= 0) goto 600
2100    else
2101      call MPI_Reduce(griduncn, 0, grid_size3d, mp_sp, MPI_SUM, id_root, &
2102           & mp_comm_used, mp_ierr)
2103    end if
2104
2105    if ((WETDEP).and.(ldirect.gt.0)) then
2106      call MPI_Reduce(wetgriduncn, wetgriduncn0, grid_size2d, mp_cp, MPI_SUM, id_root, &
2107           & mp_comm_used, mp_ierr)
2108      if (mp_ierr /= 0) goto 600
2109    end if
2110
2111    if ((DRYDEP).and.(ldirect.gt.0)) then
2112      call MPI_Reduce(drygriduncn, drygriduncn0, grid_size2d, mp_cp, MPI_SUM, id_root, &
2113           & mp_comm_used, mp_ierr)
2114      if (mp_ierr /= 0) goto 600
2115    end if
2116
2117    if (mp_measure_time) call mpif_mtime('commtime',1)
2118
2119    goto 601
2120
2121600 write(*,*) "mpi_mod> mp_ierr \= 0", mp_ierr
2122    stop
2123
2124601 end subroutine mpif_tm_reduce_grid_nest
2125
2126
2127  subroutine mpif_mtime(ident,imode)
2128!***********************************************************************
2129! Measure CPU/Wall time in various parts of the code
2130!
2131! VARIABLES
2132!   ident        character, identifies routine to measure
2133!   imode        integer, 0:start clock(s)  1: stop clock(s)
2134!
2135!***********************************************************************
2136    implicit none
2137
2138    character(LEN=*), intent(in) :: ident
2139    integer, intent(in) :: imode
2140
2141!***********************************************************************
2142
2143    select case(ident)
2144
2145    case ('timemanager')
2146      if (imode.eq.0) then
2147        call cpu_time(tm_tot_beg)
2148      else
2149        call cpu_time(tm_tot_end)
2150        tm_tot_total = tm_tot_total + (tm_tot_end-tm_tot_beg)
2151      end if
2152
2153    case ('wetdepo')
2154      if (imode.eq.0) then
2155        mp_wetdepo_wtime_beg = mpi_wtime()
2156        call cpu_time(mp_wetdepo_time_beg)
2157      else
2158        mp_wetdepo_wtime_end = mpi_wtime()
2159        call cpu_time(mp_wetdepo_time_end)
2160
2161        mp_wetdepo_wtime_total = mp_wetdepo_wtime_total + &
2162             &(mp_wetdepo_wtime_end - mp_wetdepo_wtime_beg)
2163        mp_wetdepo_time_total = mp_wetdepo_time_total + &
2164             &(mp_wetdepo_time_end - mp_wetdepo_time_beg)
2165      end if
2166
2167    case ('advance')
2168      if (imode.eq.0) then
2169        mp_advance_wtime_beg = mpi_wtime()
2170      else
2171        mp_advance_wtime_end = mpi_wtime()
2172
2173        mp_advance_wtime_total = mp_advance_wtime_total + &
2174             &(mp_advance_wtime_end - mp_advance_wtime_beg)
2175      end if
2176
2177    case ('getfields')
2178      if (imode.eq.0) then
2179        mp_getfields_wtime_beg = mpi_wtime()
2180        call cpu_time(mp_getfields_time_beg)
2181      else
2182        mp_getfields_wtime_end = mpi_wtime()
2183        call cpu_time(mp_getfields_time_end)
2184
2185        mp_getfields_wtime_total = mp_getfields_wtime_total + &
2186             &(mp_getfields_wtime_end - mp_getfields_wtime_beg)
2187        mp_getfields_time_total = mp_getfields_time_total + &
2188             &(mp_getfields_time_end - mp_getfields_time_beg)
2189      end if
2190
2191    case ('partloop1')
2192      if (imode.eq.0) then
2193        call cpu_time(tm_nploop_beg)
2194      else
2195        call cpu_time(tm_nploop_end)
2196        tm_nploop_total = tm_nploop_total + (tm_nploop_end - tm_nploop_beg)
2197      end if
2198
2199    case ('conccalc')
2200      if (imode.eq.0) then
2201        call cpu_time(mp_conccalc_time_beg)
2202      else
2203        call cpu_time(mp_conccalc_time_end)
2204        mp_conccalc_time_total = mp_conccalc_time_total + mp_conccalc_time_end - &
2205             &mp_conccalc_time_beg
2206      end if
2207
2208    case ('rootonly')
2209      if (imode.eq.0) then
2210        call cpu_time(mp_root_time_beg)
2211        mp_root_wtime_beg = mpi_wtime()
2212      else
2213        call cpu_time(mp_root_time_end)
2214        mp_root_wtime_end = mpi_wtime()
2215        mp_root_time_total = mp_root_time_total + &
2216             &(mp_root_time_end - mp_root_time_beg)
2217        mp_root_wtime_total = mp_root_wtime_total + &
2218             &(mp_root_wtime_end - mp_root_wtime_beg)
2219      end if
2220
2221    case ('iotime')
2222      if (imode.eq.0) then
2223        mp_io_wtime_beg = mpi_wtime()
2224        call cpu_time(mp_io_time_beg)
2225      else
2226        mp_io_wtime_end = mpi_wtime()
2227        call cpu_time(mp_io_time_end)
2228
2229        mp_io_wtime_total = mp_io_wtime_total + (mp_io_wtime_end - &
2230             & mp_io_wtime_beg)
2231        mp_io_time_total = mp_io_time_total + (mp_io_time_end - &
2232             & mp_io_time_beg)
2233      end if
2234
2235    case ('verttransform')
2236      if (imode.eq.0) then
2237        mp_vt_wtime_beg = mpi_wtime()
2238        call cpu_time(mp_vt_time_beg)
2239      else
2240        mp_vt_wtime_end = mpi_wtime()
2241        call cpu_time(mp_vt_time_end)
2242
2243        mp_vt_wtime_total = mp_vt_wtime_total + (mp_vt_wtime_end - &
2244             & mp_vt_wtime_beg)
2245        mp_vt_time_total = mp_vt_time_total + (mp_vt_time_end - &
2246             & mp_vt_time_beg)
2247      end if
2248
2249    case ('readwind')
2250      if (imode.eq.0) then
2251        call cpu_time(mp_readwind_time_beg)
2252        mp_readwind_wtime_beg = mpi_wtime()
2253      else
2254        call cpu_time(mp_readwind_time_end)
2255        mp_readwind_wtime_end = mpi_wtime()
2256
2257        mp_readwind_time_total = mp_readwind_time_total + &
2258             &(mp_readwind_time_end - mp_readwind_time_beg)
2259        mp_readwind_wtime_total = mp_readwind_wtime_total + &
2260             &(mp_readwind_wtime_end - mp_readwind_wtime_beg)
2261      end if
2262
2263    case ('commtime')
2264      if (imode.eq.0) then
2265        call cpu_time(mp_comm_time_beg)
2266        mp_comm_wtime_beg = mpi_wtime()
2267      else
2268        call cpu_time(mp_comm_time_end)
2269        mp_comm_wtime_end = mpi_wtime()
2270        mp_comm_time_total = mp_comm_time_total + &
2271             &(mp_comm_time_end - mp_comm_time_beg)
2272        mp_comm_wtime_total = mp_comm_wtime_total + &
2273             &(mp_comm_wtime_end - mp_comm_wtime_beg)
2274      end if
2275
2276    case ('flexpart')
2277      if (imode.eq.0) then
2278        mp_total_wtime_beg=mpi_wtime()
2279      else
2280        mp_total_wtime_end=mpi_wtime()
2281        mp_total_wtime_total = mp_total_wtime_end-mp_total_wtime_beg
2282      end if
2283
2284    case default
2285      write(*,*) 'mpi_mod::mpif_mtime> WARNING: unknown case identifier'
2286
2287    end select
2288
2289
2290  end subroutine mpif_mtime
2291
2292
2293  subroutine mpif_finalize
2294!***********************************************************************
2295! Finalize MPI
2296! Optionally print detailed time diagnostics
2297!***********************************************************************
2298    implicit none
2299
2300    integer :: ip !,j,r
2301
2302!***********************************************************************
2303
2304    IF (mp_measure_time) THEN
2305      do ip=0, mp_np-1
2306        call MPI_BARRIER(MPI_COMM_WORLD, mp_ierr)
2307
2308        if (ip == mp_pid) then
2309          write(*,FMT='(72("#"))')
2310          write(*,FMT='(A60,I3)') 'STATISTICS FOR MPI PROCESS:', mp_pid
2311          if (ip == 0) then
2312            write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR ROOT-ONLY (PID 0) &
2313                 &SECTIONS: ', mp_root_time_total
2314            write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR ROOT-ONLY (PID 0) &
2315                 &SECTIONS:', mp_root_wtime_total
2316          end if
2317          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR FLEXPART: ' &
2318               &, mp_total_wtime_total
2319          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR MPI COMMUNICATIONS:',&
2320               & mp_comm_time_total
2321          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR MPI COMMUNICATIONS:',&
2322               & mp_comm_wtime_total
2323          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR MPI BARRIERS:',&
2324               & mp_barrier_time_total
2325          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR MPI BARRIERS:',&
2326               & mp_barrier_wtime_total
2327          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR TIMEMANAGER:',&
2328               & tm_tot_total
2329          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR TIMEMANAGER/NUMPART LOOP:',&
2330               & tm_nploop_total
2331          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR ADVANCE:',&
2332               & mp_advance_wtime_total
2333          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR GETFIELDS:',&
2334               & mp_getfields_wtime_total
2335          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR GETFIELDS:',&
2336               & mp_getfields_time_total
2337          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR READWIND:',&
2338               & mp_readwind_wtime_total
2339          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR READWIND:',&
2340               & mp_readwind_time_total
2341          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR FILE IO:',&
2342               & mp_io_wtime_total
2343          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR FILE IO:',&
2344               & mp_io_time_total
2345          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR WETDEPO:',&
2346               & mp_wetdepo_wtime_total
2347          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR WETDEPO:',&
2348               & mp_wetdepo_time_total
2349          write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR CONCCALC:',&
2350               & mp_conccalc_time_total
2351          ! write(*,FMT='(A60,TR1,F9.2)') 'TOTAL WALL TIME FOR VERTTRANSFORM:',&
2352          !      & mp_vt_wtime_total
2353          ! write(*,FMT='(A60,TR1,F9.2)') 'TOTAL CPU TIME FOR VERTTRANSFORM:',&
2354          !      & mp_vt_time_total
2355! NB: the 'flush' function is possibly a gfortran-specific extension
2356          call flush()
2357        end if
2358      end do
2359    end if
2360
2361! This call to barrier is for correctly formatting output
2362    call MPI_BARRIER(MPI_COMM_WORLD, mp_ierr)
2363
2364    if (lroot.and.mp_measure_time) then
2365      write(*,FMT='(72("#"))')
2366      WRITE(*,*) "To turn off output of time measurements, set "
2367      WRITE(*,*) "    mp_measure_time=.false."
2368      WRITE(*,*) "in file mpi_mod.f90"
2369      write(*,FMT='(72("#"))')
2370    end if
2371
2372! j=mp_pid*nvar_async
2373! In the implementation with 3 fields, the processes may have posted
2374! MPI_Irecv requests that should be cancelled here
2375! if (.not.lmp_sync) then
2376!   r=mp_pid*nvar_async
2377!   do j=r,r+nvar_async-1
2378!     call MPI_Cancel(j,mp_ierr)
2379!     if (mp_ierr /= 0) write(*,*) '#### mpif_finalize::MPI_Cancel> ERROR ####'
2380!   end do
2381! end if
2382
2383    call MPI_FINALIZE(mp_ierr)
2384    if (mp_ierr /= 0) then
2385      write(*,*) '#### mpif_finalize::MPI_FINALIZE> MPI ERROR ', mp_ierr, ' ####'
2386      stop
2387    end if
2388
2389
2390    end subroutine mpif_finalize
2391
2392
2393    subroutine get_lun(my_lun)
2394!***********************************************************************
2395! get_lun:
2396!   Starting from 100, get next free logical unit number
2397!***********************************************************************
2398
2399      implicit none
2400
2401      integer, intent(inout) :: my_lun
2402      integer, save :: free_lun=100
2403      logical :: exists, iopen
2404
2405!***********************************************************************
2406
2407      loop1: do
2408        inquire(UNIT=free_lun, EXIST=exists, OPENED=iopen)
2409        if (exists .and. .not.iopen) exit loop1
2410        free_lun = free_lun+1
2411      end do loop1
2412      my_lun = free_lun
2413
2414    end subroutine get_lun
2415
2416
2417    subroutine write_data_dbg(array_in, array_name, tstep, ident)
2418!***********************************************************************
2419! Write one-dimensional arrays to file (for debugging purposes)
2420!***********************************************************************
2421      implicit none
2422
2423      real, intent(in), dimension(:) :: array_in
2424      integer, intent(in) :: tstep
2425      integer :: lios
2426      character(LEN=*), intent(in) :: ident, array_name
2427
2428      character(LEN=8) :: c_ts
2429      character(LEN=40) :: fn_1, fn_2
2430
2431!***********************************************************************
2432
2433      write(c_ts, FMT='(I8.8,BZ)') tstep
2434      fn_1='-'//trim(adjustl(c_ts))//'-'//trim(ident)
2435      write(c_ts, FMT='(I2.2,BZ)') mp_np
2436      fn_2= trim(adjustl(array_name))//trim(adjustl(fn_1))//'-np'//trim(adjustl(c_ts))//'.dat'
2437
2438      call get_lun(dat_lun)
2439      open(UNIT=dat_lun, FILE=fn_2, IOSTAT=lios, ACTION='WRITE', &
2440           FORM='UNFORMATTED', STATUS='REPLACE')
2441      write(UNIT=dat_lun, IOSTAT=lios) array_in
2442      close(UNIT=dat_lun)
2443
2444    end subroutine write_data_dbg
2445
2446
2447  subroutine set_fields_synthetic()
2448!*******************************************************************************
2449! DESCRIPTION
2450!   Set all meteorological fields to synthetic (usually constant/homogeneous)
2451!   values.
2452!   Used for validation and error-checking
2453!
2454! NOTE
2455!   This version uses asynchronious communications.
2456!
2457! VARIABLES
2458!   
2459!
2460!
2461!*******************************************************************************
2462    use com_mod
2463
2464    implicit none
2465
2466    integer :: li=1, ui=2 ! wfmem indices (i.e, operate on entire field)
2467
2468    if (.not.lmp_sync) ui=3
2469
2470
2471! Variables transferred at initialization only
2472!*********************************************
2473!    readclouds=readclouds_
2474    oro(:,:)=0.0
2475    excessoro(:,:)=0.0
2476    lsm(:,:)=0.0
2477    xlanduse(:,:,:)=0.0
2478!    wftime
2479!    numbwf
2480!    nmixz
2481!    height
2482
2483! Time-varying fields:
2484    uu(:,:,:,li:ui) = 10.0
2485    vv(:,:,:,li:ui) = 0.0
2486    uupol(:,:,:,li:ui) = 10.0
2487    vvpol(:,:,:,li:ui)=0.0
2488    ww(:,:,:,li:ui)=0.
2489    tt(:,:,:,li:ui)=300.
2490    rho(:,:,:,li:ui)=1.3
2491    drhodz(:,:,:,li:ui)=.0
2492    tth(:,:,:,li:ui)=0.0
2493    qvh(:,:,:,li:ui)=1.0
2494    qv(:,:,:,li:ui)=1.0
2495
2496    pv(:,:,:,li:ui)=1.0
2497    clouds(:,:,:,li:ui)=0
2498
2499    clwc(:,:,:,li:ui)=0.0
2500    ciwc(:,:,:,li:ui)=0.0
2501 
2502! 2D fields
2503
2504    cloudsh(:,:,li:ui)=0
2505    vdep(:,:,:,li:ui)=0.0
2506    ps(:,:,:,li:ui)=1.0e5
2507    sd(:,:,:,li:ui)=0.0
2508    tcc(:,:,:,li:ui)=0.0
2509    tt2(:,:,:,li:ui)=300.
2510    td2(:,:,:,li:ui)=250.
2511    lsprec(:,:,:,li:ui)=0.0
2512    convprec(:,:,:,li:ui)=0.0
2513    ustar(:,:,:,li:ui)=1.0
2514    wstar(:,:,:,li:ui)=1.0
2515    hmix(:,:,:,li:ui)=10000.
2516    tropopause(:,:,:,li:ui)=10000.
2517    oli(:,:,:,li:ui)=0.01
2518 
2519  end subroutine set_fields_synthetic
2520
2521end module mpi_mod
Note: See TracBrowser for help on using the repository browser.
hosted by ZAMG