source: CPL/oasis3-mct/branches/OASIS3-MCT_2.0_branch/examples/tutorial/model2.F90_TP @ 4775

Last change on this file since 4775 was 4775, checked in by aclsce, 4 years ago
  • Imported oasis3-mct from Cerfacs svn server (not suppotred anymore).

The version has been extracted from https://oasis3mct.cerfacs.fr/svn/branches/OASIS3-MCT_2.0_branch/oasis3-mct@1818

File size: 13.1 KB
Line 
1!------------------------------------------------------------------------
2! Copyright 2010, CERFACS, Toulouse, France.
3! All rights reserved. Use is subject to OASIS3 license terms.
4!=============================================================================
5!
6PROGRAM model2
7  !
8  ! Use for netCDF library
9  USE netcdf
10  ! Use for OASIS communication library
11  USE mod_oasis
12  !
13  IMPLICIT NONE
14
15  INCLUDE 'mpif.h'
16
17  !
18  ! By default OASIS3 exchanges data in double precision.
19  ! To exchange data in single precision with OASIS3,
20  ! the coupler has to be compiled with CPP key "use_realtype_single"
21  ! and the model with CPP key "NO_USE_DOUBLE_PRECISION"
22#ifdef NO_USE_DOUBLE_PRECISION
23  INTEGER, PARAMETER :: wp = SELECTED_REAL_KIND(6,37)   ! real
24#elif USE_DOUBLE_PRECISION
25  INTEGER, PARAMETER :: wp = SELECTED_REAL_KIND(12,307) ! double
26#endif
27  !
28  CHARACTER(len=30), PARAMETER :: data_filename='grid_model2.nc'
29  ! Component name (6 characters) same as in the namcouple
30  CHARACTER(len=6)   :: comp_name = 'model2'
31  CHARACTER(len=128) :: comp_out ! name of the output log file
32  CHARACTER(len=3)   :: chout
33  !
34  ! Global grid parameters :
35  INTEGER :: nlon, nlat     ! dimensions in the 2 directions of space
36  INTEGER :: ntot           ! total dimension
37  INTEGER :: il_paral_size
38  INTEGER :: nc ! number of corners in the (i,j) plan
39  INTEGER :: indi_beg, indi_end, indj_beg, indj_end
40  !
41  DOUBLE PRECISION, DIMENSION(:,:), POINTER   :: globalgrid_lon,globalgrid_lat
42  DOUBLE PRECISION, DIMENSION(:,:,:), POINTER :: globalgrid_clo,globalgrid_cla
43  DOUBLE PRECISION, DIMENSION(:,:), POINTER   :: globalgrid_srf
44  INTEGER, DIMENSION(:,:), POINTER            :: indice_mask ! mask, 0 == valid point, 1 == masked point
45  !
46  INTEGER :: mype, npes ! rank and number of pe
47  INTEGER :: localComm  ! local MPI communicator and Initialized
48  INTEGER :: comp_id    ! component identification
49  !
50  INTEGER, DIMENSION(:), ALLOCATABLE :: il_paral ! Decomposition for each proc
51  !
52  INTEGER :: ierror, rank, w_unit
53  INTEGER :: i, j
54  !
55  ! Names of exchanged Fields
56  CHARACTER(len=8), PARAMETER :: var_name1 = 'FRECVATM' ! 8 characters field received by the atmosphere from the ocean
57  CHARACTER(len=8), PARAMETER :: var_name2 = 'FSENDATM' ! 8 characters field sent by the atmosphere to the ocean
58  !
59  ! Used in oasis_def_var and oasis_def_var
60  INTEGER                   :: var_id(2)
61  INTEGER                   :: var_nodims(2)
62  INTEGER                   :: var_type
63  !
64  REAL (kind=wp), PARAMETER :: field_ini = -1. ! initialisation of received fields
65  !
66  INTEGER               ::  ib
67  INTEGER, PARAMETER    ::  il_nb_time_steps = 12 ! number of time steps
68  INTEGER, PARAMETER    ::  delta_t = 1800       ! time step
69  !
70  INTEGER                :: il_flag          ! Flag for grid writing
71  !
72  INTEGER                :: itap_sec ! Time
73  !
74  ! Grid parameter definition
75  INTEGER                :: part_id  ! use to connect the partition to the variables
76                                     ! in oasis_def_var
77  INTEGER                :: var_actual_shape(4) ! local dimensions of the arrays to the pe
78                                                ! 2 x field rank (= 4 because fields are of rank = 2)
79  !
80  ! Exchanged local fields arrays
81  ! used in routines oasis_put and oasis_get
82  REAL (kind=wp), POINTER :: field1_recv(:,:)
83  REAL (kind=wp), POINTER :: field2_send(:,:)
84  !
85  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
86  !   INITIALISATION
87  !++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
88  !
89  !!!!!!!!!!!!!!!!! OASIS_INIT !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
90  !
91  CALL oasis_init_comp (comp_id, comp_name, ierror )
92  IF (ierror /= 0) THEN
93      WRITE(0,*) 'oasis_init_comp abort by model2 compid ',comp_id
94      CALL oasis_abort(comp_id,comp_name,'Problem at line 98')
95  ENDIF
96  !
97  ! Unit for output messages : one file for each process
98  CALL MPI_Comm_Rank ( MPI_COMM_WORLD, rank, ierror )
99  IF (ierror /= 0) THEN
100      WRITE(0,*) 'MPI_Comm_Rank abort by model2 compid ',comp_id
101      CALL oasis_abort(comp_id,comp_name,'Problem at line 105')
102  ENDIF
103  !
104  w_unit = 100 + rank
105  WRITE(chout,'(I3)') w_unit
106  comp_out=comp_name//'.out_'//chout
107  !
108  OPEN(w_unit,file=TRIM(comp_out),form='formatted')
109  WRITE (w_unit,*) '-----------------------------------------------------------'
110  WRITE (w_unit,*) 'MPI_COMM_WORLD is :',MPI_COMM_WORLD
111  WRITE (w_unit,*) TRIM(comp_name), ' Running with reals compiled as kind =',wp
112  WRITE (w_unit,*) 'I am component ', TRIM(comp_name), ' rank :',rank
113  WRITE (w_unit,*) '----------------------------------------------------------'
114  CALL flush(w_unit)
115  ! 
116  !!!!!!!!!!!!!!!!! OASIS_GET_LOCALCOMM !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
117  !
118  CALL oasis_get_localcomm ( localComm, ierror )
119  IF (ierror /= 0) THEN
120      WRITE (w_unit,*) 'oasis_get_localcomm abort by model2 compid ',comp_id
121      CALL oasis_abort(comp_id,comp_name,'Problem at line 125')
122  ENDIF
123  !
124  ! Get MPI size and rank
125  CALL MPI_Comm_Size ( localComm, npes, ierror )
126  IF (ierror /= 0) THEN
127      WRITE(w_unit,*) 'MPI_comm_size abort by model2 compid ',comp_id
128      CALL oasis_abort(comp_id,comp_name,'Problem at line 132')
129  ENDIF
130  !
131  CALL MPI_Comm_Rank ( localComm, mype, ierror )
132  IF (ierror /= 0) THEN
133      WRITE (w_unit,*) 'MPI_Comm_Rank abort by model2 compid ',comp_id
134      CALL oasis_abort(comp_id,comp_name,'Problem at line 138')
135  ENDIF
136  !
137  WRITE(w_unit,*) 'I am the ', TRIM(comp_name), ' ', 'comp', comp_id, 'local rank', mype
138  WRITE (w_unit,*) 'Number of processors :',npes
139  WRITE (w_unit,*) 'Local communicator :', localComm
140  CALL flush(w_unit)
141  !
142  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
143  !  GRID DEFINITION
144  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
145  !
146  ! Reading netcdf file with pre-defined variable names
147  !
148  ! Reading dimensions of the grid
149  CALL read_dimgrid(nlon,nlat,data_filename,w_unit)
150  nc=4
151  !
152  ! Allocation
153  ALLOCATE(globalgrid_lon(nlon,nlat), STAT=ierror )
154  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating globalgrid_lon'
155  ALLOCATE(globalgrid_lat(nlon,nlat), STAT=ierror )
156  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating globalgrid_lat'
157  ALLOCATE(globalgrid_clo(nlon,nlat,nc), STAT=ierror )
158  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating globalgrid_clo'
159  ALLOCATE(globalgrid_cla(nlon,nlat,nc), STAT=ierror )
160  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating globalgrid_cla'
161  ALLOCATE(globalgrid_srf(nlon,nlat), STAT=ierror )
162  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating globalgrid_srf'
163  ALLOCATE(indice_mask(nlon,nlat), STAT=ierror )
164  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating indice_mask'
165  !
166  ! Reading of the longitudes, latitudes, longitude and latitudes of the corners, mask of the grid
167  CALL read_grid(nlon,nlat,nc,data_filename,w_unit, &
168                 globalgrid_lon,globalgrid_lat, &
169                 globalgrid_clo,globalgrid_cla, &
170                 globalgrid_srf, &
171                 indice_mask)
172  !
173  ! (Global) grid definition for OASIS3
174  ! Writing of the file grids.nc and masks.nc by the processor 0 from the grid read in
175  !
176  IF (mype == 0) THEN
177      !
178      ! Mask inversion to follow (historical) OASIS3 convention (0=not masked;1=masked)
179      WHERE(indice_mask == 1)
180          indice_mask=0
181      ELSEWHERE
182          indice_mask=1
183      END WHERE
184      !
185      CALL oasis_start_grids_writing(il_flag)
186      CALL oasis_write_grid('lmdz', nlon, nlat, globalgrid_lon, globalgrid_lat)
187      CALL oasis_write_corner('lmdz', nlon, nlat, 4, globalgrid_clo, globalgrid_cla)
188      call oasis_write_area('lmdz', nlon, nlat, globalgrid_srf)
189      CALL oasis_write_mask('lmdz', nlon, nlat, indice_mask(:,:))
190      CALL oasis_terminate_grids_writing()
191  ENDIF
192  WRITE(w_unit,*) 'After grids writing'
193  call flush(w_unit)
194  !
195  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
196  !  PARTITION DEFINITION
197  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ !
198  !
199  ! Definition of the partition of the grid (calling oasis_def_partition)
200  ntot=nlon*nlat
201#ifdef DECOMP_APPLE
202  il_paral_size = 3
203#elif defined DECOMP_BOX
204  il_paral_size = 5
205#endif
206  ALLOCATE(il_paral(il_paral_size))
207  WRITE(w_unit,*) 'After allocate il_paral, il_paral_size', il_paral_size
208  call flush(w_unit)
209  !
210  CALL decomp_def (il_paral,il_paral_size,nlon,nlat,mype,npes,w_unit)
211  WRITE(w_unit,*) 'After decomp_def, il_paral = ', il_paral(:)
212  call flush(w_unit)
213  ! The data are exchanged in the global grid so you do not need to pass
214  ! isize to oasis_def_partition
215  CALL oasis_def_partition (part_id, il_paral, ierror)
216  !
217  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
218  ! DEFINITION OF THE LOCAL FIELDS 
219  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
220  !
221  !!!!!!!!!!!!!!!!!! OASIS_DEF_VAR !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
222  !
223  ! Define transient variables
224  !
225  var_nodims(1) = 2    ! Rank of the field array is 2
226  var_nodims(2) = 1    ! Bundles always 1 for OASIS3
227  var_type = OASIS_Real
228  !
229  var_actual_shape(1) = 1
230  var_actual_shape(2) = il_paral(3)
231  var_actual_shape(3) = 1
232#ifdef DECOMP_APPLE
233  var_actual_shape(4) = 1
234#elif defined DECOMP_BOX
235  var_actual_shape(4) = il_paral(4)
236#endif
237  ! Declaration of the field associated with the partition of the grid
238  CALL oasis_def_var (var_id(1),var_name1, part_id, &
239     var_nodims, OASIS_In, var_actual_shape, var_type, ierror)
240  IF (ierror /= 0) THEN
241      WRITE (w_unit,*) 'oasis_def_var abort by model2 compid ',comp_id
242      CALL oasis_abort(comp_id,comp_name,'Problem at line 242')
243  ENDIF
244  !
245  CALL oasis_def_var (var_id(2),var_name2, part_id, &
246     var_nodims, OASIS_Out, var_actual_shape, var_type, ierror)
247  IF (ierror /= 0) THEN
248      WRITE (w_unit,*) 'oasis_def_var abort by model2 compid ',comp_id
249      CALL oasis_abort(comp_id,comp_name,'Problem at line 249')
250  ENDIF
251  !
252  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
253  !         TERMINATION OF DEFINITION PHASE
254  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
255  !  All processes involved in the coupling must call oasis_enddef;
256  !  here all processes are involved in coupling
257  !
258  !!!!!!!!!!!!!!!!!! OASIS_ENDDEF !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
259  !
260  CALL oasis_enddef ( ierror )
261  IF (ierror /= 0) THEN
262      WRITE (w_unit,*) 'oasis_enddef abort by model2 compid ',comp_id
263      CALL oasis_abort(comp_id,comp_name,'Problem at line 263')
264  ENDIF
265  !
266  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
267  ! SEND AND RECEIVE ARRAYS
268  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
269  !
270  ! Allocate the fields send and received by the model
271  !
272  ALLOCATE(field1_recv(var_actual_shape(2), var_actual_shape(4)), STAT=ierror )
273  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating field1_recv'
274  !
275  ALLOCATE(field2_send(var_actual_shape(2), var_actual_shape(4)),STAT=ierror )
276  IF ( ierror /= 0 ) WRITE(w_unit,*) 'Error allocating field2_send'
277  !
278  DEALLOCATE(il_paral)
279  !
280  !!!!!!!!!!!!!!!!!!!!!!!!!!!!OASIS_PUT/OASIS_GET !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
281  !
282  indi_beg=1 ; indi_end=nlon
283  indj_beg=((nlat/npes)*mype)+1
284  !
285  IF (mype .LT. npes - 1) THEN
286      indj_end = (nlat/npes)*(mype+1)
287  ELSE
288      indj_end = nlat
289  ENDIF
290  !
291  ! Data exchange
292  !
293  ! Time loop
294  DO ib=1, il_nb_time_steps
295    itap_sec = delta_t * (ib-1) ! Time
296    !
297    ! Get the field FRECVATM
298    field1_recv=field_ini
299    CALL oasis_get(var_id(1),itap_sec, field1_recv, ierror)
300    write(w_unit,*) 'tcx recvf1 ',itap_sec,minval(field1_recv),maxval(field1_recv)
301    IF ( ierror .NE. OASIS_Ok .AND. ierror .LT. OASIS_Recvd) THEN
302        WRITE (w_unit,*) 'oasis_get abort by model2 compid ',comp_id
303        CALL oasis_abort(comp_id,comp_name,'Problem at line 306')
304    ENDIF
305    !
306    ! Send the field FSENDATM
307    !
308    CALL function_sent(var_actual_shape(2), var_actual_shape(4), &
309                       RESHAPE(globalgrid_lon(indi_beg:indi_end,indj_beg:indj_end),&
310                               (/ var_actual_shape(2), var_actual_shape(4) /)), &
311                       RESHAPE(globalgrid_lat(indi_beg:indi_end,indj_beg:indj_end),&
312                               (/ var_actual_shape(2), var_actual_shape(4) /)), &
313                                field2_send,ib)
314
315    !                 
316    write(w_unit,*) 'tcx sendf2 ',itap_sec,minval(field2_send),maxval(field2_send)
317    CALL oasis_put(var_id(2),itap_sec, field2_send, ierror)
318    IF ( ierror .NE. OASIS_Ok .AND. ierror .LT. OASIS_Sent) THEN
319        WRITE (w_unit,*) 'oasis_put abort by model2 compid ',comp_id
320        CALL oasis_abort(comp_id,comp_name,'Problem at line 317')
321    ENDIF
322    !
323  ENDDO
324  !
325  WRITE (w_unit,*) 'End of the program'
326  CALL flush(w_unit)
327  CLOSE(w_unit)
328  !
329  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
330  !         TERMINATION
331  !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
332  !
333  !!!!!!!!!!!!!!!!!! OASIS_TERMINATE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
334  !
335  ! Collective call to terminate the coupling exchanges
336  !
337  CALL oasis_terminate (ierror)
338  IF (ierror /= 0) THEN
339      WRITE (w_unit,*) 'oasis_terminate abort by model2 compid ',comp_id
340      CALL oasis_abort(comp_id,comp_name,'Problem at line 337')
341  ENDIF
342  !
343END PROGRAM MODEL2
344!
Note: See TracBrowser for help on using the repository browser.