/[MITgcm]/MITgcm_contrib/nesting_sannino/nest_driver/driver_nesting.F
ViewVC logotype

Contents of /MITgcm_contrib/nesting_sannino/nest_driver/driver_nesting.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.2 - (show annotations) (download)
Sun Nov 28 03:26:49 2010 UTC (15 years ago) by jmc
Branch: MAIN
CVS Tags: HEAD
Changes since 1.1: +0 -0 lines
FILE REMOVED
plit driver_nesting.F into "main.F" and "interpolation_p2c.F".

1 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
2 C
3 PROGRAM NEST_DRIVER
4 C !DESCRIPTION:
5 C Routine that manages the MPI communication between the CHILD
6 C and PARENT models. It performs also the necessary
7 C interpolations from PARENT2CHILD and CHILD2PARENT.
8 C
9 C ver 1.0 by G. Sannino, V. Ruggiero, A. Carillo, P. Heimbach
10 C
11 C First application described in:
12 C Sannino G.,Herrmann, Carillo, Rupolo, Ruggiero, Artale, Heimbach, 2009:
13 C An eddy-permitting model of the Mediterranean Sea with a two-way grid
14 C refinement at Gibraltar.
15 C Ocean Modelling, 30(1), 56-72, doi: 10.1016/j.ocemod.2009.06.002
16
17 C
18
19 C !LOCAL INPUT VARIABLES:
20 c ---------------------------------------------------------------------------------
21 c NST_LEV_TOT :: Total nesting levels (1 for only one nesting)
22 c NST_LEV :: Number of the actual nesting
23 c NCPUs_CHLD :: Number of CPUs used for the CHILD model
24 c at NST_LEV nesting level
25 c NCPUs_PRNT :: Number of CPUs used for the PARENT model
26 c at NST_LEV nesting level
27 c nSxC,nSyF :: Domain decomposition used for CHILD
28 c nSxP,nSyP :: Domain decomposition used for PARENT
29 c OLX,OLY :: Domain dec. overlapping (same for both models)
30 c NrC,NyC,NxC :: Dimension PARENT model
31 c NrF,NyF,NxF :: Dimension CHILD model
32 c WesterB :: Western (i) PARENT index where start the refinement
33 c EasterB :: Eastern (i) PARENT index where finish the refinement
34 c dirNEST :: Directory where are stored the geometry data of both models
35 c ---------------------------------------------------------------------------------
36 CEOP
37 implicit none
38 c--------------------------------------------------------
39 c INPUT VARIABLE DEFINITION
40 c--------------------------------------------------------
41 INTEGER :: NST_LEV_TOT, NST_LEV, NCPUs_PRNT
42 INTEGER :: Count_Lev
43 PARAMETER (NST_LEV_TOT = 1) !Number of Total Nesting Levels
44 PARAMETER (NST_LEV = 1) !Which level am I?
45
46 INTEGER :: NCPUs_CHLD(NST_LEV_TOT)
47 INTEGER :: MSTR_DRV(NST_LEV_TOT)
48 INTEGER :: MSTR_PRNT(NST_LEV_TOT)
49 INTEGER :: MSTR_CHLD(NST_LEV_TOT)
50
51 PARAMETER (NCPUs_PRNT = 40)
52
53 DATA NCPUs_CHLD / 20 /
54 c--------------------------------------------------------
55 INTEGER :: nSxC,nSyC !Domain decomposition CHILD
56 INTEGER :: nSxP,nSyP !Domain decomposition PARENT
57 PARAMETER (nSxC = 4 , nSyC = 5)
58 PARAMETER (nSxP = 8 , nSyP = 5)
59 c--------------------------------------------------------
60 INTEGER :: OLY,OLX !Domain decomposition overlapping
61 c !(the same for both models)
62 PARAMETER (OLX = 3, OLY = 3)
63 c--------------------------------------------------------
64 INTEGER :: NrP,NxP,NyP
65 INTEGER :: NrC,NxC,NyC
66 INTEGER :: IM_C,JM_C
67 INTEGER :: IM_P,JM_P
68 INTEGER :: IndI,IndJ
69 INTEGER :: IndI_P(nSxP*nSyP),IndJ_P(nSxP*nSyP)
70 INTEGER :: IndI_C(nSxC*nSyC),IndJ_C(nSxC*nSyC)
71
72 INTEGER :: WesternB,EasternB
73 c--------------------------------------------------------
74 PARAMETER (NrP=42, NyP=120,NxP = 400) !PARENT model
75 PARAMETER (NrC=42, NyC=105,NxC = 140) !CHILD model
76 c--------------------------------------------------------
77 PARAMETER (WesternB = 43,EasternB=90)
78 c--------------------------------------------------------
79 CHARACTER :: dirNEST*80
80 c--------------------------------------------------------
81 PARAMETER (dirNEST ="/home/sannino/NESTING/")
82 c--------------------------------------------------------
83 INCLUDE 'mpif.h'
84 INTEGER :: ierr,rank,size,npd
85 INTEGER :: irank,isize
86 INTEGER :: color
87 INTEGER :: istatus,NEST_comm
88 INTEGER :: from,whm,status(MPI_STATUS_SIZE),st_count
89 INTEGER :: I,J,K,II,JJ,Irec,III,JJJ,KK,ICONT
90 INTEGER :: iVar,Indx,Jndx
91 INTEGER :: J1,J2,JJ1,JJ2
92 INTEGER :: I_START,I_END,I_STEP
93 REAL*4 :: XF,YF,XP1,YP1,XP2,YP2,YP3
94 REAL*8 :: TRANSPORT_WEST,TRANSPORT_EAST
95 CHARACTER*10 :: c2i(30)
96 c----------------------------------------------------
97 c Define PARENT Model Geometry
98 c----------------------------------------------------
99 c REAL*4 Xu_P(NxP,NyP)
100 REAL*4 :: Yu_P(NxP,NyP)
101 REAL*4 :: Xv_P(NxP,NyP)
102 REAL*4 :: Yv_P(NxP,NyP)
103 REAL*4 :: Xo_P(NxP,NyP)
104 REAL*4 :: Yo_P(NxP,NyP)
105 REAL*4 :: Xg_P(NxP,NyP)
106 REAL*4 :: Yg_P(NxP,NyP)
107 REAL*4 :: hFacW_P(NxP,NyP,NrP)
108 REAL*4 :: hFacS_P(NxP,NyP,NrP)
109 REAL*4 :: RAC_P(NxP,NyP)
110 REAL*4 :: RAW_P(NxP,NyP)
111 REAL*4 :: RAS_P(NxP,NyP)
112 REAL*4 :: hFacC_P(NxP,NyP,NrP)
113 REAL*4 :: DEEP_P(NxP,NyP,NrP)
114 REAL*4 :: INV_VOL_C_P(NxP,NyP,NrP)
115 REAL*4 :: INV_VOL_S_P(NxP,NyP,NrP)
116 REAL*4 :: INV_VOL_W_P(NxP,NyP,NrP)
117 c----------------------------------------------------
118 c Define CHILD Model Geometry
119 c----------------------------------------------------
120 REAL*4 :: Xu_C(NxC,NyC)
121 REAL*4 :: Yu_C(NxC,NyC)
122 REAL*4 :: Xv_C(NxC,NyC)
123 REAL*4 :: Yv_C(NxC,NyC)
124 REAL*4 :: Xo_C(NxC,NyC)
125 REAL*4 :: Yo_C(NxC,NyC)
126 REAL*4 :: Xg_C(NxC,NyC)
127 REAL*4 :: Yg_C(NxC,NyC)
128 REAL*4 :: hFacW_C(NxC,NyC,NrC)
129 REAL*4 :: hFacS_C(NxC,NyC,NrC)
130 REAL*4 :: RAC_C(NxC,NyC)
131 REAL*4 :: RAW_C(NxC,NyC)
132 REAL*4 :: RAS_C(NxC,NyC)
133 REAL*4 :: hFacC_C(NxC,NyC,NrC)
134 REAL*4 :: DEEP_C(NxC,NyC,NrC)
135 c----------------------------------------------------
136 c Define relative (PARENT-->CHILD) indicies
137 c----------------------------------------------------
138 INTEGER :: P2C_U(NyC)
139 c
140 INTEGER :: P2C_linU(NyC) !Linear interp.
141 INTEGER :: WO3_linU(NyC) !Linear interp. !Which Of 3
142 c
143 INTEGER :: P2C_linV(NyC) !Linear interp.
144 INTEGER :: WO3_linV(NyC) !Linear interp. !Which Of 3
145 c
146 INTEGER :: P2C_V(NyC) !Linear interp.
147 INTEGER :: P2C_o(NyC) !Linear interp.
148 INTEGER :: P2C1_V(NyC) !BiLinear interp.
149 INTEGER :: P2C2_V(NyC) !BiLinear interp.
150 INTEGER :: P2C1_o(NyC) !BiLinear interp.
151 INTEGER :: P2C2_o(NyC) !BiLinear interp.
152 c----------------------------------------------------
153 c Define relative (CHILD-->PARENT) indicies
154 c----------------------------------------------------
155 INTEGER I_C2P(9,NxP,NyP)
156 INTEGER J_C2P(9,NxP,NyP)
157 c----------------------------------------------------
158 c Define CHILD model variable
159 c----------------------------------------------------
160 c _____________ (1) WesternB (2) EasternB
161 c |
162 REAL*8 :: U_C1(NyC,NrC,2)
163 REAL*8 :: V_C1(NyC,NrC,2)
164 REAL*8 :: T_C1(NyC,NrC,2)
165 REAL*8 :: S_C1(NyC,NrC,2)
166 REAL*8 :: ETA_C1(NyC,NrC,2)
167 INTEGER :: MSIZE
168 c
169 REAL*8 :: U_C2(NyC,NrC,2)
170 REAL*8 :: V_C2(NyC,NrC,2)
171 REAL*8 :: T_C2(NyC,NrC,2)
172 REAL*8 :: S_C2(NyC,NrC,2)
173 REAL*8 :: ETA_C2(NyC,NrC,2)
174 c
175 REAL*8,allocatable :: VAR_C1(:,:,:,:)
176 c
177 REAL*8 :: DIFF_U(NyC,NrC,2)
178 REAL*8 :: DIFF_V(NyC,NrC,2)
179 REAL*8 :: DIFF_T(NyC,NrC,2)
180 REAL*8 :: DIFF_S(NyC,NrC,2)
181 REAL*8 :: DIFF_ETA(NyC,NrC,2)
182 c----------------------------------------------------
183 c Define PARENT model variable
184 c----------------------------------------------------
185 REAL*8 :: VAR3D_P(NxP,NyP,NrP,15)
186 REAL*8 :: VAR2D_P(NxP,NyP,4)
187
188 INTEGER :: ONOFF
189 INTEGER :: index_var3D,index_var2D
190 c---------------------------------------------------------------|
191 c (1) U || (2) V || (3) T || (4) S |
192 c---------------------------------------------------------------|
193 c (5) gU || (6) gV || (7) gT || (8) gS |
194 c---------------------------------------------------------------|
195 c (9) gUNm1 || (10) gVNm1 || (11) gTNm1 || (12) gSNm1 |
196 c---------------------------------------------------------------|
197 c (13) totPhiHyd || (14) IVDConvCount || |
198 c---------------------------------------------------------------|
199 c (15) etaN || (16) etaH || (17) phiHydLow |
200 c---------------------------------------------------------------|
201 c (18) etaNm1 || (19) etaHm1|| |
202 c---------------------------------------------------------------|
203 c---------------------------------------------------------------|
204 c Define Global Variables to Exchange |
205 c---------------------------------------------------------------|
206 REAL*8,allocatable :: globalPA (:,:,:,:) !(6,NyP,NrP,5)
207 REAL*8 :: globalP1(6,NyP,NrP)
208 REAL*8 :: globalP2(6,NyP,NrP)
209 REAL*8 :: globalP3(6,NyP,NrP)
210 REAL*8 :: globalP4(6,NyP,NrP)
211 REAL*8 :: globalP5(6,NyP,NrP)
212 REAL*8 :: globalP6(6,NyP,NrP)
213 REAL*8 :: globalP7(6,NyP,NrP)
214 REAL*8 :: globalP8(6,NyP,NrP)
215 REAL*8 :: globalP9(6,NyP,NrP)
216 REAL*8 :: globalP10(6,NyP,NrP)
217 REAL*8 :: globalP11(6,NyP,NrP)
218 REAL*8 :: globalP12(6,NyP,NrP)
219 REAL*8 :: globalP13(6,NyP,NrP)
220 REAL*8 :: globalP14(6,NyP,NrP)
221 c
222 INTEGER :: index
223 c----------------------------------------------------
224 c Define Global Variables to Exchange
225 c----------------------------------------------------
226 REAL*8 :: globalC3D(NxC,NyC,NrC,15)
227 c |___________ 15 fields
228 c
229
230 REAL*8 globalC2D(NxC,NyC,4)
231 c |___________ 4 fields
232 c
233 c
234 REAL*8,allocatable :: globalC3D_a(:,:,:,:),globalC2D_a(:,:,:)
235 c
236 INTEGER :: indexF,index1F
237 REAL*4 :: AREA_VOL
238 INTEGER :: vstart,vstop,VCONT,VCONTP(0:3)
239 c----------------------------------------------------
240 c MPI starts here
241 c----------------------------------------------------
242 call MPI_Init(ierr)
243 call MPI_Comm_size(MPI_COMM_WORLD,size,ierr)
244 call MPI_Comm_rank(MPI_COMM_WORLD,rank,ierr)
245 npd=size-(NCPUs_PRNT+NCPUs_CHLD(1))
246 if(rank.lt.npd) color=0
247 call MPI_COMM_SPLIT (MPI_COMM_WORLD, color,0,
248 & NEST_COMM,ierr)
249 c
250 call MPI_Comm_size(NEST_COMM,isize,ierr)
251 call MPI_Comm_rank(NEST_COMM,irank,ierr)
252 c--------------------------------------------------------
253 c COMPUTE MASTER VALUES
254 c--------------------------------------------------------
255 MSTR_DRV(1) = 0
256 c
257 MSTR_PRNT(1) = npd
258 MSTR_CHLD(1) = NCPUs_PRNT + npd
259 c
260 DO Count_Lev = 2, NST_LEV_TOT
261 MSTR_DRV(Count_Lev) = MSTR_CHLD(Count_Lev-1) +
262 & NCPUs_CHLD(Count_Lev - 1)
263 c
264 MSTR_CHLD(Count_Lev) = MSTR_DRV(Count_Lev) + 1
265 MSTR_PRNT(Count_Lev) = MSTR_CHLD(Count_Lev-1)
266 ENDDO
267 c
268 vstart = 1+rank*(nSxP/MSTR_PRNT(1))
269 vstop = (1+rank)*(nSxP/MSTR_PRNT(1))
270 VCONT = (nSxP/npd)*(nSyP)*rank-1
271 VCONTP(rank) = VCONT
272 c--------------------------------------------------------
273 c COMPUTE DOMAIN DECOMPOSITION
274 c--------------------------------------------------------
275 if(rank.eq.0) then
276 c
277 IM_C = int(NxC/nSxC)
278 JM_C = int(NyC/nSyC)
279 c
280 IM_P = int(NxP/nSxP)
281 JM_P = int(NyP/nSyP)
282 c
283 indexF = IM_C*JM_C*NrC*15
284 index1F = IM_C*JM_C*4
285 c
286 indexF = (JM_C+OLY+OLY)*NrC*2*5
287 c
288 index_var3D = IM_P*JM_P*NrP
289 index_var2D = IM_P*JM_P
290 c
291 ICONT = 0
292 DO I = 1,nSxP
293 DO J = 1,nSyP
294 ICONT = ICONT + 1
295 IndI_P(ICONT) = IM_P*(I-1)
296 IndJ_P(ICONT) = JM_P*(J-1)
297 END DO
298 END DO
299 c
300 ICONT = 0
301 DO I = 1,nSxC
302 DO J = 1,nSyC
303 ICONT = ICONT + 1
304 IndI_C(ICONT) = IM_C*(I-1)
305 IndJ_C(ICONT) = JM_C*(J-1)
306 END DO
307 END DO
308 c
309 ALLOCATE (globalPA (6,JM_P,NrP,5))
310 index=6*JM_P*NrP*5
311 c
312 ALLOCATE(globalC3D_a(IM_C,JM_C,NrC,15))
313 ALLOCATE(globalC2D_a(IM_C,JM_C,4))
314 c
315 ALLOCATE (VAR_C1((JM_C+OLY+OLY),NrC,2,5))
316 c--------------------------------------------------------
317 c WARNING
318 c--------------------------------------------------------
319 print*,'*************************************'
320 print*,' have you update geometric files?'
321 print*,' in ./CHILD e ./PARENT'
322 print*,'*************************************'
323 c--------------------------------------------------------
324 c PARENT MODEL
325 c--------------------------------------------------------
326 print*,' [1] Read PARENT model geometry'
327 c----------------------------------------------------
328 c XC & YC
329 c----------------------------------------------------
330 MSIZE = NxP*NyP
331 c
332 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
333 & convert='big_endian',
334 & file=trim(dirNEST)//'/PARENT/XC.data',
335 & form='unformatted')
336 c
337 read (1,REC=1) Xo_P(:,:)
338 close(1)
339 c
340 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
341 & convert='big_endian',
342 & file=trim(dirNEST)//'/PARENT/YC.data',
343 & form='unformatted')
344 c
345 read (1,REC=1) Yo_P(:,:)
346 close(1)
347 c----------------------------------------------------
348 c XG & YG
349 c----------------------------------------------------
350 MSIZE = NxP*NyP
351 c
352 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
353 & convert='big_endian',
354 & file=trim(dirNEST)//'/PARENT/XG.data',
355 & form='unformatted')
356
357 read (1,REC=1) Xg_P(:,:)
358 close(1)
359 c
360 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
361 & convert='big_endian',
362 & file=trim(dirNEST)//'/PARENT/YG.data',
363 & form='unformatted')
364 c
365 read (1,REC=1) Yg_P(:,:)
366 close(1)
367 c----------------------------------------------------
368 c Yu
369 c----------------------------------------------------
370 DO J = 1,NyP
371 DO I = 1,NxP
372 c Xu_P(I,J) = Xg_P(I,J)
373 Yu_P(I,J) = Yo_P(I,J)
374 ENDDO
375 ENDDO
376 c----------------------------------------------------
377 c Xv & Yv
378 c----------------------------------------------------
379 DO J = 1,NyP
380 DO I = 1,NxP
381 Xv_P(I,J) = Xo_P(I,J)
382 Yv_P(I,J) = Yg_P(I,J)
383 ENDDO
384 ENDDO
385 c----------------------------------------------------
386 c hFacC
387 c----------------------------------------------------
388 MSIZE = NxP*NyP*NrP
389 c
390 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
391 & convert='big_endian',
392 & file=trim(dirNEST)//'/PARENT/hFacC.data',
393 & form='unformatted')
394
395 read (1,REC=1) hFacC_P(:,:,:)
396 close(1)
397 c----------------------------------------------------
398 c hFacW
399 c----------------------------------------------------
400 MSIZE = NxP*NyP*NrP
401 c
402 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
403 & convert='big_endian',
404 & file=trim(dirNEST)//'/PARENT/hFacW.data',
405 & form='unformatted')
406
407 read (1,REC=1) hFacW_P(:,:,:)
408 close(1)
409 c----------------------------------------------------
410 c hFacS
411 c----------------------------------------------------
412 MSIZE = NxP*NyP*NrP
413 c
414 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
415 & convert='big_endian',
416 & file=trim(dirNEST)//'/PARENT/hFacS.data',
417 & form='unformatted')
418
419 read (1,REC=1) hFacS_P(:,:,:)
420 close(1)
421 c----------------------------------------------------
422 c RAC
423 c----------------------------------------------------
424 MSIZE = NxP*NyP
425 c
426 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
427 & convert='big_endian',
428 & file=trim(dirNEST)//'/PARENT/RAC.data',
429 & form='unformatted')
430
431 read (1,REC=1) RAC_P(:,:)
432 close(1)
433 c----------------------------------------------------
434 c RAW
435 c----------------------------------------------------
436 MSIZE = NxP*NyP
437 c
438 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
439 & convert='big_endian',
440 & file=trim(dirNEST)//'/PARENT/RAW.data',
441 & form='unformatted')
442
443 read (1,REC=1) RAW_P(:,:)
444 close(1)
445 c----------------------------------------------------
446 c RAS
447 c----------------------------------------------------
448 MSIZE = NxP*NyP
449 c
450 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
451 & convert='big_endian',
452 & file=trim(dirNEST)//'/PARENT/RAS.data',
453 & form='unformatted')
454
455 read (1,REC=1) RAS_P(:,:)
456 close(1)
457 c----------------------------------------------------
458 c MASK x PARENT
459 c----------------------------------------------------
460 DO K = 1,NrP
461 DO J = 1,NyP
462 DO I = 1,NxP
463 DEEP_P(i,j,k) = 0.
464 IF (hFacC_P(i,j,k).ne.0) then
465 DEEP_P(I,J,K) = 1.
466 ENDIF
467 ENDDO
468 ENDDO
469 ENDDO
470 c----------------------------------------------------
471 c 1/Volume (C)
472 c----------------------------------------------------
473 DO K = 1,NrP
474 DO J = 1,NyP
475 DO I = 1,NxP
476 INV_VOL_C_P(I,J,K) = 1.
477 IF ((RAC_P(I,J)*hFacC_P(I,J,K)).ne.0.) THEN
478 INV_VOL_C_P(I,J,K) = 1./(RAC_P(I,J)*hFacC_P(I,J,K))
479 ENDIF
480 ENDDO
481 ENDDO
482 ENDDO
483 c----------------------------------------------------
484 c 1/Volume (W)
485 c----------------------------------------------------
486 DO K = 1,NrP
487 DO J = 1,NyP
488 DO I = 1,NxP
489 INV_VOL_W_P(I,J,K) = 1.
490 IF ((RAW_P(I,J)*hFacW_P(I,J,K)).ne.0.) THEN
491 INV_VOL_W_P(I,J,K) = 1./(RAW_P(I,J)*hFacW_P(I,J,K))
492 ENDIF
493 ENDDO
494 ENDDO
495 ENDDO
496 c----------------------------------------------------
497 c 1/Volume (S)
498 c----------------------------------------------------
499 DO K = 1,NrP
500 DO J = 1,NyP
501 DO I = 1,NxP
502 INV_VOL_S_P(I,J,K) = 1.
503 IF ((RAS_P(I,J)*hFacS_P(I,J,K)).ne.0.) THEN
504 INV_VOL_S_P(I,J,K) = 1./(RAS_P(I,J)*hFacS_P(I,J,K))
505 ENDIF
506 ENDDO
507 ENDDO
508 ENDDO
509 c--------------------------------------------------------
510 c CHILD MODEL
511 c--------------------------------------------------------
512 print*,' [2] Read CHILD model geometry'
513 c----------------------------------------------------
514 c XC & YC
515 c----------------------------------------------------
516 MSIZE = NxC*NyC
517 c
518 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
519 & convert='big_endian',
520 & file=trim(dirNEST)//'/CHILD/XC.data',
521 & form='unformatted')
522
523 read (1,REC=1) Xo_C(:,:)
524 close(1)
525 c
526 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
527 & convert='big_endian',
528 & file=trim(dirNEST)//'/CHILD/YC.data',
529 & form='unformatted')
530 c
531 read (1,REC=1) Yo_C(:,:)
532 close(1)
533 c----------------------------------------------------
534 c XG & YG
535 c----------------------------------------------------
536 MSIZE = NxC*NyC
537 c
538 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
539 & convert='big_endian',
540 & file=trim(dirNEST)//'/CHILD/XG.data',
541 & form='unformatted')
542
543 read (1,REC=1) Xg_C(:,:)
544 close(1)
545 c
546 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
547 & convert='big_endian',
548 & file=trim(dirNEST)//'/CHILD/YG.data',
549 & form='unformatted')
550
551 read (1,REC=1) Yg_C(:,:)
552 close(1)
553 c----------------------------------------------------
554 c Xu & Yu
555 c----------------------------------------------------
556 DO J = 1,NyC
557 DO I = 1,NxC
558 Xu_C(I,J) = Xg_C(I,J)
559 Yu_C(I,J) = Yo_C(I,J)
560 ENDDO
561 ENDDO
562 c----------------------------------------------------
563 c Xv & Yv
564 c----------------------------------------------------
565 DO J = 1,NyC
566 DO I = 1,NxC
567 Xv_C(I,J) = Xo_C(I,J)
568 Yv_C(I,J) = Yg_C(I,J)
569 ENDDO
570 ENDDO
571 c----------------------------------------------------
572 c hFacC
573 c----------------------------------------------------
574 MSIZE = NxC*NyC*NrC
575 c
576 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
577 & convert='big_endian',
578 & file=trim(dirNEST)//'/CHILD/hFacC.data',
579 & form='unformatted')
580
581 read (1,REC=1) hFacC_C(:,:,:)
582 close(1)
583 c----------------------------------------------------
584 c hFacW
585 c----------------------------------------------------
586 MSIZE = NxC*NyC*NrC
587 c
588 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
589 & convert='big_endian',
590 & file=trim(dirNEST)//'/CHILD/hFacW.data',
591 & form='unformatted')
592
593 read (1,REC=1) hFacW_C(:,:,:)
594 close(1)
595 c----------------------------------------------------
596 c hFacC
597 c----------------------------------------------------
598 MSIZE = NxC*NyC*NrC
599 c
600 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
601 & convert='big_endian',
602 & file=trim(dirNEST)//'/CHILD/hFacS.data',
603 & form='unformatted')
604
605 read (1,REC=1) hFacS_C(:,:,:)
606 close(1)
607 c----------------------------------------------------
608 c RAC
609 c----------------------------------------------------
610 MSIZE = NxC*NyC
611 c
612 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
613 & convert='big_endian',
614 & file=trim(dirNEST)//'/CHILD/RAC.data',
615 & form='unformatted')
616
617 read (1,REC=1) RAC_C(:,:)
618 close(1)
619 c----------------------------------------------------
620 c RAW
621 c----------------------------------------------------
622 MSIZE = NxC*NyC
623 c
624 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
625 & convert='big_endian',
626 & file=trim(dirNEST)//'/CHILD/RAW.data',
627 & form='unformatted')
628
629 read (1,REC=1) RAW_C(:,:)
630 close(1)
631 c----------------------------------------------------
632 c RAS
633 c----------------------------------------------------
634 MSIZE = NxC*NyC
635 c
636 open(unit=1,ACCESS='direct',RECL=MSIZE,STATUS='OLD',
637 & convert='big_endian',
638 & file=trim(dirNEST)//'/CHILD/RAS.data',
639 & form='unformatted')
640
641 read (1,REC=1) RAS_C(:,:)
642 close(1)
643 c----------------------------------------------------
644 c MASK x CHILD
645 c----------------------------------------------------
646 DO K = 1,NrC
647 DO J = 1,NyC
648 DO I = 1,NxC
649 DEEP_C(i,j,k) = 0.
650 IF (hFacC_C(i,j,k).ne.0) then
651 DEEP_C(I,J,K) = 1.
652 ENDIF
653 ENDDO
654 ENDDO
655 ENDDO
656 c----------------------------------------------------
657 c
658 c __/________ ___________
659 c | | | | ||
660 c > o > | | |
661 c |__/__|_____|
662 c | | |
663 c > o > |
664 c |_____|_____|_____|_____|
665 c
666 c
667 c----------------------------------------------------
668 print*,' [3] Compute J index P-->C'
669 c--------------------------------------------------------
670 c Compute J indicies for mapping P->C (I)
671 c--------------------------------------------------------
672 I = 1
673 II = WesternB
674 c
675 DO J = 1,NyC
676 P2C_U(J) = 0.
677 DO JJ = 1,NyP-1
678 YF = Yg_C(I,J)
679 YP1 = Yg_P(II,JJ)
680 YP2 = Yg_P(II,JJ+1)
681 IF (YF.ge.YP1.and.YF.lt.YP2) THEN
682 P2C_U(J) = JJ
683 ENDIF
684 ENDDO
685 ENDDO
686 c--------------------------------------------------------
687 c Compute J indicies for mapping P->C (II)
688 c--------------------------------------------------------
689 I = 1
690 II = WesternB
691 c
692 DO J = 1,NyC
693 P2C_linU(J) = 0.
694 DO JJ = 1,NyP-1
695 YF = Yu_C(I,J)
696 YP1 = Yu_P(II,JJ)
697 YP2 = Yu_P(II,JJ+1)
698 IF (YF.ge.YP1.and.YF.lt.YP2) THEN
699 P2C_linU(J) = JJ
700 ENDIF
701 ENDDO
702 ENDDO
703 c--------------------------------------------------------
704 c Compute J indicies for mapping P->C (III)
705 c--------------------------------------------------------
706 I = 1
707 II = WesternB
708 c
709 DO J = 1,NyC
710 DO JJ = 1,NyP-1
711 YF = Yu_C(I,J)
712 YP1 = Yu_P(II,JJ)
713 IF (YF.eq.YP1) THEN
714 WO3_linU(J) = 0
715 if (J+1.le.NyC) WO3_linU(J+1) = 1
716 if (J+2.le.NyC) WO3_linU(J+2) = 2
717 ENDIF
718 ENDDO
719 ENDDO
720 c--------------------Lower bound
721 DO J = 1,NyC
722 DO JJ = 1,NyP-1
723 YF = Yu_C(I,J)
724 YP1 = Yu_P(II,JJ)
725 IF (YF.eq.YP1) THEN
726 WO3_linU(J) = 0
727 if (J-1.gt.0) WO3_linU(J-1) = 2
728 if (J-2.gt.0) WO3_linU(J-2) = 1
729 GOTO 2345
730 ENDIF
731 ENDDO
732 ENDDO
733 2345 CONTINUE
734 c--------------------Upper bound
735 DO J = NyC,1,-1
736 DO JJ = 1,NyP-1
737 YF = Yu_C(I,J)
738 YP1 = Yu_P(II,JJ)
739 IF (YF.eq.YP1) THEN
740 WO3_linU(J) = 0
741 if (J+1.le.NyC) WO3_linU(J+1) = 1
742 if (J+2.le.NyC) WO3_linU(J+2) = 2
743 GOTO 2346
744 ENDIF
745 ENDDO
746 ENDDO
747 2346 CONTINUE
748 c--------------------------------------------------------
749 c Compute J indicies for mapping P->C (IV)
750 c--------------------------------------------------------
751 I = 1
752 II = WesternB
753 c
754 DO J = 1,NyC
755 P2C_linV(J) = 0.
756 DO JJ = 1,NyP-1
757 YF = Yv_C(I,J)
758 YP1 = Yv_P(II,JJ)
759 YP2 = Yv_P(II,JJ+1)
760 IF (YF.ge.YP1.and.YF.lt.YP2) THEN
761 P2C_linV(J) = JJ
762 ENDIF
763 ENDDO
764 ENDDO
765 c--------------------------------------------------------
766 c Compute J indicies for mapping P->C (V)
767 c--------------------------------------------------------
768 I = 1
769 II = WesternB
770 c
771 DO J = 1,NyC
772 DO JJ = 1,NyP-1
773 YF = Yv_C(I,J)
774 YP1 = Yv_P(II,JJ)
775 IF (YF.eq.YP1) THEN
776 WO3_linV(J) = 0
777 if (J+1.le.NyC) WO3_linV(J+1) = 1
778 if (J+2.le.NyC) WO3_linV(J+2) = 2
779 ENDIF
780 ENDDO
781 ENDDO
782 c--------------------Lower bound
783 DO J = 1,NyC
784 DO JJ = 1,NyP-1
785 YF = Yv_C(I,J)
786 YP1 = Yv_P(II,JJ)
787 IF (YF.eq.YP1) THEN
788 WO3_linV(J) = 0
789 if (J-1.gt.0) WO3_linV(J-1) = 2
790 if (J-2.gt.0) WO3_linV(J-2) = 1
791 GOTO 23451
792 ENDIF
793 ENDDO
794 ENDDO
795 23451 CONTINUE
796 c--------------------Upper bound
797 DO J = NyC,1,-1
798 DO JJ = 1,NyP-1
799 YF = Yv_C(I,J)
800 YP1 = Yv_P(II,JJ)
801 IF (YF.eq.YP1) THEN
802 WO3_linV(J) = 0
803 if (J+1.le.NyC) WO3_linV(J+1) = 1
804 if (J+2.le.NyC) WO3_linV(J+2) = 2
805 GOTO 23461
806 ENDIF
807 ENDDO
808 ENDDO
809 23461 CONTINUE
810 c--------------------------------------------------------
811 c Compute J indicies for mapping P->C (V)
812 c--------------------------------------------------------
813 print*,' [5] Compute J index P-->C for (o)'
814 I = 1
815 II = WesternB
816 c
817 DO J = 1,NyC
818 P2C_o(J) = 0.
819 DO JJ = 1,NyP-1
820 YF = Yo_C(I,J)
821 YP1 = Yg_P(II,JJ)
822 YP2 = Yg_P(II,JJ+1)
823 IF (YF.gt.YP1.and.YF.lt.YP2) THEN
824 P2C_o(J) = JJ
825 ENDIF
826 ENDDO
827 ENDDO
828 c--------------------------------------------------------
829 c Compute J indicies for mapping P->C (VI)
830 c--------------------------------------------------------
831 print*,' [6] Compute J index P-->C for (v bilinear)'
832 I = 1
833 II = WesternB
834 c
835 DO J = 1,NyC
836 DO JJ = 2,NyP-1
837 YF = Yv_C(I,J)
838 YP1 = Yv_P(II,JJ)
839 YP2 = Yv_P(II,JJ+1)
840 YP3 = Yv_P(II,JJ-1)
841 c
842 IF (YF.ge.YP1.and.YF.lt.YP2) THEN
843 P2C1_V(J) = JJ
844 P2C2_V(J) = JJ+1
845 ENDIF
846 ENDDO
847 ENDDO
848 c--------------------------------------------------------
849 c Look for the 9 CHILD indicies in PARENT grid cell
850 c--------------------------------------------------------
851 print*,' [8] Compute I J index C-->P for (o)'
852 c
853 DO J = 1,NyP
854 DO I = 1,NxP
855 I_C2P(:,I,J) = 0
856 J_C2P(:,I,J) = 0
857 c
858 DO JJ = 1,NyC
859 DO II = 1,NxC
860 IF (Xo_C(II,JJ).eq.Xo_P(I,J).and.
861 & Yo_C(II,JJ).eq.Yo_P(I,J)) then
862
863 KK = 0
864 DO JJJ = JJ-1,JJ+1
865 DO III = II-1,II+1
866 KK = kk +1
867 if (III.lt.1.or.III.gt.NxC) cycle
868 if (JJJ.lt.1.or.JJJ.gt.NyC) cycle
869 I_C2P(KK,I,J) = III
870 J_C2P(KK,I,J) = JJJ
871 ENDDO
872 ENDDO
873 ENDIF
874
875 ENDDO
876 ENDDO
877 c
878 ENDDO
879 ENDDO
880 ENDIF
881 c--------------------------------------------------------
882 c Broadcast all the above variables
883 c--------------------------------------------------------
884 call MPI_BCAST(I_C2P,9*NxP*NyP,MPI_INTEGER,
885 & 0,NEST_COMM,ierr)
886 call MPI_BCAST(J_C2P,9*NxP*NyP,MPI_INTEGER,
887 & 0,NEST_COMM,ierr)
888
889 call MPI_BCAST(RAC_C,NxC*NyC,MPI_REAL,
890 & 0,NEST_COMM,ierr)
891 call MPI_BCAST(hFacC_C,NxC*NyC*NrC,MPI_REAL,
892 & 0,NEST_COMM,ierr)
893 call MPI_BCAST(INV_VOL_C_P,NxP*NyP*NrP,MPI_REAL,
894 & 0,NEST_COMM,ierr)
895
896 call MPI_BCAST(RAW_C,NxC*NyC,MPI_REAL,
897 & 0,NEST_COMM,ierr)
898 call MPI_BCAST(hFacW_C,NxC*NyC*NrC,MPI_REAL,
899 & 0,NEST_COMM,ierr)
900 call MPI_BCAST(INV_VOL_W_P,NxP*NyP*NrP,MPI_REAL,
901 & 0,NEST_COMM,ierr)
902
903 call MPI_BCAST(RAS_C,NxC*NyC,MPI_REAL,
904 & 0,NEST_COMM,ierr)
905 call MPI_BCAST(hFacS_C,NxC*NyC*NrC,MPI_REAL,
906 & 0,NEST_COMM,ierr)
907 call MPI_BCAST(INV_VOL_S_P,NxP*NyP*NrP,MPI_REAL,
908 & 0,NEST_COMM,ierr)
909 c
910 call MPI_BCAST(DEEP_C,NxC*NyC*NrC,MPI_REAL,
911 & 0,NEST_COMM,ierr)
912 call MPI_BCAST(RAC_P,NxP*NyP,MPI_REAL,
913 & 0,NEST_COMM,ierr)
914 c
915 call MPI_BCAST(IM_P,1,MPI_INTEGER,
916 & 0,NEST_COMM,ierr)
917 call MPI_BCAST(JM_P,1,MPI_INTEGER,
918 & 0,NEST_COMM,ierr)
919 call MPI_BCAST(index_var3D,1,MPI_INTEGER,
920 & 0,NEST_COMM,ierr)
921 call MPI_BCAST(index_var2D,1,MPI_INTEGER,
922 & 0,NEST_COMM,ierr)
923 c
924 call MPI_BCAST(DEEP_P,NxP*NyP*NrP,MPI_REAL,
925 & 0,NEST_COMM,ierr)
926 call MPI_BCAST(hFacS_P,NxP*NyP*NrP,MPI_REAL,
927 & 0,NEST_COMM,ierr)
928 call MPI_BCAST(hFacC_P,NxP*NyP*NrP,MPI_REAL,
929 & 0,NEST_COMM,ierr)
930 call MPI_BCAST(hFacW_P,NxP*NyP*NrP,MPI_REAL,
931 & 0,NEST_COMM,ierr)
932
933 c--------------------------------------------------------
934 if(rank.eq.0) then
935 c--------------------------------------------------------
936 DO K = 1,NrP
937 DO J = 1,NyP
938 DO I = WesternB+1,EasternB-1
939 cc WesternB side
940
941 DO II = 1,9
942 IF (I_C2P(II,I,J).eq.0) cycle
943 IF (J_C2P(II,I,J).eq.0) cycle
944 c
945 Indx = I_C2P(II,I,J)
946 Jndx = J_C2P(II,I,J)
947 ENDDO
948 ENDDO
949 ENDDO
950 ENDDO
951 c---------------------------------------------------------
952 ONOFF=0
953 endif
954 c--------------------------------------------------------
955 c BEGIN MAIN LOOP
956 c--------------------------------------------------------
957 do
958 if(rank.eq.0) then
959 c--------------------------------------------------------
960 c (1) READ FROM PARENT MODEL
961 c--------------------------------------------------------
962 ICONT=1
963 DO WHILE(ICONT.le.nSxP*nSyP)
964 from= MPI_ANY_SOURCE
965
966 call MPI_RECV (globalPA, index, MPI_REAL8,
967 & FROM, 3000,
968 & MPI_COMM_World, status,ierr)
969 c
970 ICONT=ICONT+1
971 c
972 whm=status(MPI_SOURCE)-MSTR_PRNT(NST_LEV)+1
973 c
974 call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr)
975 c
976 DO II = 1,6
977 IF (globalPA(II,1,1,1).ne.-999.) THEN
978 globalP1(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
979 & globalPA(II,1:JM_P,:,1)
980 globalP2(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
981 & globalPA(II,1:JM_P,:,2)
982 globalP3(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
983 & globalPA(II,1:JM_P,:,3)
984 globalP4(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
985 & globalPA(II,1:JM_P,:,4)
986 globalP5(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
987 & globalPA(II,1:JM_P,:,5)
988 ENDIF
989
990 ENDDO
991 ENDDO
992 c--------------------------------------------------------
993 c Start interpolation for CHILD
994 c--------------------------------------------------------
995 CALL INTERPOLATION_P2C (
996 & globalP1,globalP2,globalP3,globalP4,globalP5,
997 & NxP,NyP,NrP,
998 & NxC,NyC,NrC,
999 $ WesternB,EasternB,
1000 $ P2C_U,P2C_V,P2C_o,P2C1_V,P2C2_V,P2C1_o,P2C2_o,
1001 $ P2C_linU,WO3_linU,P2C_linV,WO3_linV,
1002 $ Xv_C,Yv_C,Xv_P,Yv_P,
1003 $ T_C1,S_C1,U_C1,V_C1,ETA_C1,
1004 $ DEEP_C,DEEP_P
1005 & )
1006 c==============================================================
1007 c Open Files from PARENT MODEL
1008 c==============================================================
1009 ICONT=1
1010
1011 do while(ICONT.le.nSxP*nSyP)
1012 from= MPI_ANY_SOURCE
1013
1014 call MPI_RECV (globalPA, index, MPI_REAL8,
1015 & FROM, 3000,
1016 & MPI_COMM_World, status,ierr)
1017
1018
1019
1020 ICONT=ICONT+1
1021
1022 whm=status(MPI_SOURCE)-MSTR_PRNT(NST_LEV)+1
1023
1024 call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr)
1025
1026 DO II = 1,6
1027 IF (globalPA(II,1,1,1).ne.-999.) THEN
1028 globalP1(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
1029 & globalPA(II,1:JM_P,:,1)
1030 globalP2(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
1031 & globalPA(II,1:JM_P,:,2)
1032 globalP3(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
1033 & globalPA(II,1:JM_P,:,3)
1034 globalP4(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
1035 & globalPA(II,1:JM_P,:,4)
1036 globalP5(II,1+IndJ_P(whm):JM_P+IndJ_P(whm),:) =
1037 & globalPA(II,1:JM_P,:,5)
1038 ENDIF
1039 ENDDO
1040
1041 end do
1042 c--------------------------------------------------------
1043 c Start inteprolation for CHILD
1044 c--------------------------------------------------------
1045 CALL INTERPOLATION_P2C (
1046 & globalP1,globalP2,globalP3,globalP4,globalP5,
1047 & NxP,NyP,NrP,
1048 & NxC,NyC,NrC,
1049 $ WesternB,EasternB,
1050 $ P2C_U,P2C_V,P2C_o,P2C1_V,P2C2_V,P2C1_o,P2C2_o,
1051 $ P2C_linU,WO3_linU,P2C_linV,WO3_linV,
1052 $ Xv_C,Yv_C,Xv_P,Yv_P,
1053 $ T_C2,S_C2,U_C2,V_C2,ETA_C2,
1054 $ DEEP_C,DEEP_P
1055 & )
1056
1057
1058
1059 c==============================================================
1060 c Temporal Interpolation OBCs x CHILD MODEL
1061 c==============================================================
1062 c 0 1200
1063 c ---+--.--.--+---- Parent
1064 c
1065 c |--|--|--
1066 c 0 800
1067 c 400
1068 c------------------------------------------------------------
1069 DO I = 1,2 ! WesternB & EasternB
1070 DIFF_T(:,:,I) = (T_C2(:,:,I) - T_C1(:,:,I))/3.
1071 DIFF_S(:,:,I) = (S_C2(:,:,I) - S_C1(:,:,I))/3.
1072 DIFF_U(:,:,I) = (U_C2(:,:,I) - U_C1(:,:,I))/3.
1073 DIFF_V(:,:,I) = (V_C2(:,:,I) - V_C1(:,:,I))/3.
1074 DIFF_eta(:,:,I) = (eta_C2(:,:,I) - eta_C1(:,:,I))/3.
1075 ENDDO
1076 c-------------------------------------------------------------
1077 c
1078 c Step 0 (Rec = 1 ==> WesternB)
1079 c------- (Rec = 2 ==> EasternB)
1080
1081 DO I = 1,2 !WesternB & EasternB
1082 T_C1(:,:,I) = T_C2(:,:,I) !+ DIFF_T(:,:,I)
1083 S_C1(:,:,I) = S_C2(:,:,I) !+ DIFF_S(:,:,I)
1084 U_C1(:,:,I) = U_C2(:,:,I) !+ DIFF_U(:,:,I)
1085 V_C1(:,:,I) = V_C2(:,:,I) !+ DIFF_V(:,:,I)
1086 ETA_C1(:,:,I) = ETA_C2(:,:,I) !+ DIFF_ETA(:,:,I)
1087 ENDDO
1088
1089
1090 if(ONOFF.eq.0) then
1091 c---------------------------------------------------------------------
1092 ICONT = -1
1093 DO I = 1,nSxC
1094 DO J = 1,nSyC
1095 ICONT = ICONT + 1
1096 IndI = IM_C*(I-1)
1097 IndJ = JM_C*(J-1)
1098
1099 VAR_C1(:,:,:,:) = 0.
1100 c
1101 J1 = 1+IndJ-OLY
1102 J2 = JM_C+IndJ+OLY
1103 c
1104 JJ1 = 1
1105 JJ2 = JM_C+OLY+OLY
1106 c
1107 IF(1 +IndJ-OLY.LT.0) THEN
1108 J1 = 1
1109 JJ1 = 4
1110 ENDIF
1111 c
1112 IF(JM_C+IndJ+OLY.GT.NyC) THEN
1113 J2 = NyC
1114 JJ2 = JM_C
1115 ENDIF
1116 c
1117 VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:)
1118 VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:)
1119 VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:)
1120 VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:)
1121 VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:)
1122 c
1123 call MPI_SEND (VAR_C1, indexF, MPI_REAL8,
1124 & MSTR_CHLD(NST_LEV)+ICONT, 3000,
1125 & MPI_Comm_World,ierr)
1126 c
1127 ENDDO
1128 ENDDO
1129 c----------------------------------------------------------------------
1130 cc write(*,*) 'VIC: MANDO SEGNALE DI OK AL CHILD PER INIZIALIZZARE'
1131 ONOFF=1
1132 ENDIF
1133 cc write(*,*) 'VIC: MANDO SEGNALE DI OK AL CHILD PER IL PASSO 1'
1134 c-----------------------------------------------------------------------
1135 ICONT = -1
1136 DO I = 1,nSxC
1137 DO J = 1,nSyC
1138 ICONT = ICONT + 1
1139 IndI = IM_C*(I-1)
1140 IndJ = JM_C*(J-1)
1141 c
1142 VAR_C1(:,:,:,:) = 0.
1143 c
1144 J1 = 1+IndJ-OLY
1145 J2 = JM_C+IndJ+OLY
1146 c
1147 JJ1 = 1
1148 JJ2 = JM_C+OLY+OLY
1149 c
1150 IF(1 +IndJ-OLY.LT.0) THEN
1151 J1 = 1
1152 JJ1 = 4
1153 ENDIF
1154 c
1155 IF(JM_C+IndJ+OLY.GT.NyC) THEN
1156 J2 = NyC
1157 JJ2 = JM_C
1158 ENDIF
1159 c
1160 VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:)
1161 VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:)
1162 VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:)
1163 VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:)
1164 VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:)
1165 c
1166 call MPI_SEND (VAR_C1, indexF, MPI_REAL8,
1167 & MSTR_CHLD(NST_LEV)+ICONT, 3000,
1168 & MPI_Comm_World,ierr)
1169
1170 ENDDO
1171 ENDDO
1172 c---------------------------------------------------------------------
1173 goto 8888
1174
1175 c
1176 c Step 1 (Rec = 3 ==> WesternB)
1177 c------- (Rec = 4 ==> EasternB)
1178
1179 DO I = 1,2 !WesternB & EasternB
1180 T_C1(:,:,I) = T_C2(:,:,I) !+ DIFF_T(:,:,I)
1181 S_C1(:,:,I) = S_C2(:,:,I) !+ DIFF_S(:,:,I)
1182 U_C1(:,:,I) = U_C2(:,:,I) !+ DIFF_U(:,:,I)
1183 V_C1(:,:,I) = V_C2(:,:,I) !+ DIFF_V(:,:,I)
1184 ETA_C1(:,:,I) = ETA_C2(:,:,I) !+ DIFF_ETA(:,:,I)
1185 ENDDO
1186 c----------------------------------------------------------
1187 ICONT = -1
1188 DO I = 1,nSxC
1189 DO J = 1,nSyC
1190 ICONT = ICONT + 1
1191 IndI = IM_C*(I-1)
1192 IndJ = JM_C*(J-1)
1193
1194 VAR_C1(:,:,:,:) = 0.
1195 c
1196 J1 = 1+IndJ-OLY
1197 J2 = JM_C+IndJ+OLY
1198 c
1199 JJ1 = 1
1200 JJ2 = JM_C+OLY+OLY
1201 c
1202 IF(1 +IndJ-OLY.LT.0) THEN
1203 J1 = 1
1204 JJ1 = 4
1205 ENDIF
1206 c
1207 IF(JM_C+IndJ+OLY.GT.NyC) THEN
1208 J2 = NyC
1209 JJ2 = JM_C
1210 ENDIF
1211 c
1212 VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:)
1213 VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:)
1214 VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:)
1215 VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:)
1216 VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:)
1217
1218 call MPI_SEND (VAR_C1, indexF, MPI_REAL8,
1219 & MSTR_CHLD(NST_LEV)+ICONT, 3000,
1220 & MPI_Comm_World,ierr)
1221
1222
1223 ENDDO
1224 ENDDO
1225 c-----------------------------------------------------------
1226 c
1227 c Step 2 (Rec = 5 ==> WesternB)
1228 c------- (Rec = 6 ==> EasternB)
1229
1230 DO I = 1,2 !WesternB & EasternB
1231 T_C1(:,:,I) = T_C2(:,:,I) !+ DIFF_T(:,:,I)
1232 S_C1(:,:,I) = S_C2(:,:,I) !+ DIFF_S(:,:,I)
1233 U_C1(:,:,I) = U_C2(:,:,I) !+ DIFF_U(:,:,I)
1234 V_C1(:,:,I) = V_C2(:,:,I) !+ DIFF_V(:,:,I)
1235 ETA_C1(:,:,I) = ETA_C2(:,:,I) !+ DIFF_ETA(:,:,I)
1236 ENDDO
1237 c----------------------------------------------------------
1238 ICONT = -1
1239 DO I = 1,nSxC
1240 DO J = 1,nSyC
1241 ICONT = ICONT + 1
1242 IndI = IM_C*(I-1)
1243 IndJ = JM_C*(J-1)
1244
1245 VAR_C1(:,:,:,:) = 0.
1246 c
1247 J1 = 1+IndJ-OLY
1248 J2 = JM_C+IndJ+OLY
1249 c
1250 JJ1 = 1
1251 JJ2 = JM_C+OLY+OLY
1252 c
1253 IF(1 +IndJ-OLY.LT.0) THEN
1254 J1 = 1
1255 JJ1 = 4
1256 ENDIF
1257 c
1258 IF(JM_C+IndJ+OLY.GT.NyC) THEN
1259 J2 = NyC
1260 JJ2 = JM_C
1261 ENDIF
1262 c
1263 VAR_C1(JJ1:JJ2,:,:,1) = U_C1(J1:J2,:,:)
1264 VAR_C1(JJ1:JJ2,:,:,2) = V_C1(J1:J2,:,:)
1265 VAR_C1(JJ1:JJ2,:,:,3) = T_C1(J1:J2,:,:)
1266 VAR_C1(JJ1:JJ2,:,:,4) = S_C1(J1:J2,:,:)
1267 VAR_C1(JJ1:JJ2,:,:,5) = ETA_C1(J1:J2,:,:)
1268
1269 call MPI_SEND (VAR_C1, indexF, MPI_REAL8,
1270 & MSTR_CHLD(NST_LEV)+ICONT, 3000,
1271 & MPI_Comm_World,ierr)
1272
1273
1274 ENDDO
1275 ENDDO
1276 c---------------------------------------------------------------
1277 8888 CONTINUE
1278 c--------------------------------------------------------
1279 c Close OBCs Files x CHILD MODEL
1280 c--------------------------------------------------------
1281 c----------------------------------------------------
1282 c------------- MANDO SEGNALE DI OK AL CHILD
1283 c----------------------------------------------------
1284 c--------------------------------------------------------
1285 c (1) READ FROM CHILD MODEL
1286 c--------------------------------------------------------
1287 call MPI_RECV (TRANSPORT_WEST, 1, MPI_REAL8,
1288 & MSTR_CHLD(NST_LEV), 3000,
1289 & MPI_COMM_World, status,ierr)
1290
1291
1292 call MPI_RECV (TRANSPORT_EAST, 1, MPI_REAL8,
1293 & MSTR_CHLD(NST_LEV), 3000,
1294 & MPI_COMM_World, status,ierr)
1295 c---------------------------------------------------------
1296 c---------------------------------------------------------
1297 ICONT=1
1298
1299 DO WHILE(ICONT.le.nSxC*nSyC)
1300 from= MPI_ANY_SOURCE
1301 call MPI_RECV (globalC3D_a,indexF, MPI_REAL8,
1302 & from, 3000, MPI_COMM_World, status,ierr)
1303
1304 ICONT=ICONT+1
1305
1306 whm=status(MPI_SOURCE)-MSTR_CHLD(NST_LEV)+1
1307
1308 call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr)
1309
1310 globalC3D(1+IndI_C(whm):IM_C+IndI_C(whm),
1311 & 1+IndJ_C(whm):JM_C+IndJ_C(whm),:,:)=
1312 & globalC3D_a(:,:,:,:)
1313 END DO
1314 c-----------------------------
1315 ICONT=1
1316
1317 DO WHILE(ICONT.le.nSxC*nSyC)
1318 from= MPI_ANY_SOURCE
1319 call MPI_RECV (globalC2D_a,index1F, MPI_REAL8,
1320 & from, 3000, MPI_COMM_World, status,ierr)
1321
1322 ICONT=ICONT+1
1323
1324 whm=status(MPI_SOURCE)-MSTR_CHLD(NST_LEV)+1
1325
1326 call MPI_GET_COUNT(status,MPI_REAL8,st_count,ierr)
1327
1328
1329 globalC2D(1+IndI_C(whm):IM_C+IndI_C(whm),
1330 & 1+IndJ_C(whm):JM_C+IndJ_C(whm),:)=
1331 & globalC2D_a(:,:,:)
1332 END DO
1333 ENDIF
1334
1335
1336
1337
1338 call MPI_BCAST(globalC3D,NxC*NyC*NrC*15,MPI_REAL8,
1339 & 0,NEST_COMM,ierr)
1340 call MPI_BCAST(globalC2D,NxC*NyC*4,MPI_REAL8,
1341 & 0,NEST_COMM,ierr)
1342
1343 2323 CONTINUE
1344 c=======================================================
1345 c (1) READ FROM CHILD MODEL
1346 c=======================================================
1347
1348 c=======================================================
1349 c (2) INTERPOLATIONS
1350 c=======================================================
1351
1352 c 3D VAR
1353 c--------
1354 DO iVar = 1,15 ! tipo di variabile
1355 DO K = 1,NrP
1356 DO J = 1,NyP
1357 DO I = WesternB+1,EasternB-1
1358 VAR3D_P(I,J,K,iVar) = 0. ! inizializzo
1359 c WesternB side
1360
1361 AREA_VOL = 0. !can be area or volume depend on the variable
1362
1363 SELECT CASE(iVar)
1364 CASE(1,5,9)
1365 I_START = 1
1366 I_END = 9
1367 I_STEP = 1 !3
1368 CASE(2,6,10)
1369 I_START = 1
1370 I_END = 9 !3
1371 I_STEP = 1
1372 CASE DEFAULT
1373 I_START = 1
1374 I_END = 9
1375 I_STEP = 1
1376 END SELECT
1377
1378
1379 DO II = I_START,I_END,I_STEP
1380
1381
1382 IF (I_C2P(II,I,J).eq.0) cycle
1383 IF (J_C2P(II,I,J).eq.0) cycle
1384 c
1385 Indx = I_C2P(II,I,J)
1386 Jndx = J_C2P(II,I,J)
1387 c
1388 c
1389 SELECT CASE(iVar)
1390
1391 CASE (1,5,9)
1392 VAR3D_P(I,J,K,ivar) = VAR3D_P(I,J,K,iVar) +
1393 & globalC3D(Indx,Jndx,K,iVar)*
1394 $ RAW_C(Indx,Jndx)*
1395 & hFacW_C(Indx,Jndx,K)
1396
1397 CASE (2,6,10)
1398 VAR3D_P(I,J,K,ivar) = VAR3D_P(I,J,K,iVar) +
1399 & globalC3D(Indx,Jndx,K,iVar)*
1400 $ RAS_C(Indx,Jndx)*
1401 & hFacS_C(Indx,Jndx,K)
1402
1403
1404 CASE DEFAULT
1405 VAR3D_P(I,J,K,ivar) = VAR3D_P(I,J,K,iVar) +
1406 & globalC3D(Indx,Jndx,K,iVar)*
1407 $ RAC_C(Indx,Jndx)*
1408 & hFacC_C(Indx,Jndx,K)
1409
1410 AREA_VOL = AREA_VOL +
1411 & RAC_C(Indx,Jndx)* hFacC_C(Indx,Jndx,K)
1412
1413 END SELECT
1414 ENDDO
1415 c-----------------------------------------------
1416 c Make a volume average
1417 c----------------------------------------------
1418 SELECT CASE(iVar)
1419 CASE (1,5,9)
1420 VAR3D_P(I,J,K,ivar) =
1421 & VAR3D_P(I,J,K,iVar)*
1422 & INV_VOL_W_P(I,J,K)
1423 if (hFacW_P(I,J,K).eq.0.) VAR3D_P(I,J,K,ivar)=0.
1424
1425 CASE (2,6,10)
1426 VAR3D_P(I,J,K,ivar) =
1427 & VAR3D_P(I,J,K,iVar)*
1428 & INV_VOL_S_P(I,J,K)
1429 if (hFacS_P(I,J,K).eq.0.) VAR3D_P(I,J,K,ivar)=0.
1430 CASE DEFAULT
1431 IF (AREA_VOL.ne.0.) then
1432 VAR3D_P(I,J,K,ivar) =
1433 & VAR3D_P(I,J,K,iVar)/AREA_VOL
1434 ENDIF
1435 if (hFacC_P(I,J,K).eq.0.) VAR3D_P(I,J,K,ivar)=0.
1436 END SELECT
1437 ENDDO
1438 ENDDO
1439 ENDDO
1440 ENDDO
1441
1442
1443
1444 c 2D VAR
1445 c--------
1446 DO iVar = 1,4
1447 DO J = 1,NyP
1448 DO I = WesternB+1,EasternB-1
1449 VAR2D_P(I,J,iVar) = 0.
1450 AREA_VOL = 0.
1451 DO II = 1,9
1452 IF (I_C2P(II,I,J).eq.0) cycle
1453 IF (J_C2P(II,I,J).eq.0) cycle
1454 c
1455 Indx = I_C2P(II,I,J)
1456 Jndx = J_C2P(II,I,J)
1457 c
1458 VAR2D_P(I,J,ivar) = VAR2D_P(I,J,iVar) +
1459 & globalC2D(Indx,Jndx,iVar)*
1460 $ RAC_C(Indx,Jndx)*
1461 & DEEP_C(Indx,Jndx,1)
1462
1463
1464 AREA_VOL = AREA_VOL +
1465 & RAC_C(Indx,Jndx)* DEEP_C(Indx,Jndx,1)
1466
1467 ENDDO
1468 c-----------------------------
1469 IF ((RAC_P(I,J)*DEEP_P(I,J,1)).ne.0.) then
1470 c IF (AREA_VOL.ne.0.) then
1471
1472 VAR2D_P(I,J,ivar) =
1473 & VAR2D_P(I,J,iVar)/
1474 & RAC_P(I,J)
1475 ENDIF
1476 c----------------------------
1477 ENDDO
1478 ENDDO
1479 ENDDO
1480 if(rank.eq.0) then
1481 c--------------------------------------------------------
1482 c Write Files for PARENT MODEL
1483 c--------------------------------------------------------
1484 c print*,' (*) Open Files for PARENT MODEL'
1485
1486 7575 CONTINUE
1487 c----------------------------------------------------
1488 c------------- MANDO SEGNALE DI OK AL PARENT
1489 c----------------------------------------------------
1490 call MPI_SEND (TRANSPORT_WEST, 1, MPI_REAL8,
1491 & MSTR_PRNT(NST_LEV), 3000,
1492 & MPI_Comm_World,ierr)
1493
1494 call MPI_SEND (TRANSPORT_EAST, 1, MPI_REAL8,
1495 & MSTR_PRNT(NST_LEV), 3000,
1496 & MPI_Comm_World,ierr)
1497
1498 ENDIF
1499 c---------------------------------------------------------
1500 !---------------------------------------------------------
1501 VCONT=VCONTP(rank)
1502
1503 DO I = vstart,vstop
1504 DO J = 1,nSyP
1505 VCONT = VCONT + 1
1506 IndI = IM_P*(I-1)
1507 IndJ = JM_P*(J-1)
1508 c-----------------------------------------------------------
1509 DO iVar=1,15
1510 CALL MPI_SEND (VAR3D_P(1+IndI:IM_P+IndI
1511 & ,1+IndJ:JM_P+IndJ,:,iVar),
1512 & index_var3D,MPI_REAL8,MSTR_PRNT(NST_LEV)+VCONT,
1513 & 3000,MPI_Comm_World,ierr)
1514
1515 ENDDO
1516
1517
1518
1519 DO iVar=1,4
1520 call MPI_SEND (VAR2D_P(1+IndI:IM_P+IndI
1521 & ,1+IndJ:JM_P+IndJ,iVar),
1522 & index_var2D,MPI_REAL8,MSTR_PRNT(NST_LEV)+VCONT,
1523 & 3000,MPI_Comm_World,ierr)
1524 ENDDO
1525
1526 c-----------------------------------------------------------
1527 END DO
1528 END DO
1529
1530 call MPI_BARRIER(NEST_COMM,ierr)
1531 end do
1532 c---------------------------------------------------------
1533 c=======================================================
1534 c END MAIN LOOP
1535 c=======================================================
1536 call MPI_FINALIZE(ierr)
1537 c---------------------------------------------------------
1538 !---------------------------------------------------------
1539
1540 101 FORMAT (I1)
1541
1542 STOP
1543 END
1544
1545
1546
1547
1548
1549
1550
1551
1552
1553
1554
1555
1556
1557
1558
1559 SUBROUTINE INTERPOLATION_P2C (
1560 & globalP1,globalP2,globalP3,globalP4,globalP5,
1561 & NxP,NyP,NrP,
1562 & NxC,NyC,NrC,
1563 $ WesternB,EasternB,
1564 $ P2C_U,P2C_V,P2C_o,P2C1_V,P2C2_V,P2C1_o,P2C2_o,
1565 $ P2C_linU,WO3_linU,P2C_linV,WO3_linV,
1566 $ Xv_F,Yv_F,Xv_P,Yv_P,
1567 $ T_F,S_F,U_F,V_F,ETA_F,
1568 $ DEEP_F,DEEP_P
1569 & )
1570 c----------------------------------------------------
1571 implicit none
1572 c----------------------------------------------------
1573 INTEGER :: I,J,K,II,JJ
1574 INTEGER :: WesternB,EasternB
1575 INTEGER :: NrP,NxP,NyP
1576 INTEGER :: NrC,NxC,NyC
1577 c----------------------------------------------------
1578 REAL*8 :: Fp,Fm,Fo,VEL_MEMO
1579 INTEGER :: INDC
1580 c----------------------------------------------------
1581 c Define Global Variables to Exchange
1582 c----------------------------------------------------
1583 REAL*8 :: globalP1(6,NyP,NrP)
1584 REAL*8 :: globalP2(6,NyP,NrP)
1585 REAL*8 :: globalP3(6,NyP,NrP)
1586 REAL*8 :: globalP4(6,NyP,NrP)
1587 REAL*8 :: globalP5(6,NyP,NrP)
1588 c----------------------------------------------------
1589 c Define CHILD Model Geometry
1590 c----------------------------------------------------
1591 REAL*4 :: Xu_F(NxC,NyC)
1592 REAL*4 :: Yu_F(NxC,NyC)
1593 REAL*4 :: Xv_F(NxC,NyC)
1594 REAL*4 :: Yv_F(NxC,NyC)
1595 REAL*4 :: Xo_F(NxC,NyC)
1596 REAL*4 :: Yo_F(NxC,NyC)
1597 REAL*4 :: Xg_F(NxC,NyC)
1598 REAL*4 :: Yg_F(NxC,NyC)
1599 REAL*4 :: DEEP_F(NxC,NyC,NrC)
1600 c----------------------------------------------------
1601 c Define PARENT Model Geometry
1602 c----------------------------------------------------
1603 c REAL*4 :: Xu_P(NxP,NyP)
1604 REAL*4 :: Yu_P(NxP,NyP)
1605 REAL*4 :: Xv_P(NxP,NyP)
1606 REAL*4 :: Yv_P(NxP,NyP)
1607 REAL*4 :: Xo_P(NxP,NyP)
1608 REAL*4 :: Yo_P(NxP,NyP)
1609 REAL*4 :: Xg_P(NxP,NyP)
1610 REAL*4 :: Yg_P(NxP,NyP)
1611 REAL*4 :: DEEP_P(NxP,NyP,NrP)
1612 REAL*4 :: DEPDEP
1613 c-----------------------------------------------------
1614 REAL*4 :: X1,X2,X3,X4,Y1,Y2,Y3,Y4
1615 REAL*4 :: f1,f2,f3,f4,f,x,y
1616 REAL*4 :: gammaT,gammaS,terzo,dueterzi
1617 REAL*4 :: gammaEta
1618 REAL*4 :: gammaV
1619 c----------------------------------------------------
1620 c Define INDICIES MATRIX
1621 c----------------------------------------------------
1622 INTEGER :: P2C_U(NyC) !x imposizione NETTA
1623 INTEGER :: P2C_linU(NyC) !x Lineare
1624 INTEGER :: WO3_linU(NyC) !x Lineare !Which Of 3
1625
1626 INTEGER :: P2C_linV(NyC) !x Lineare
1627 INTEGER :: WO3_linV(NyC) !x Lineare !Which Of 3
1628
1629 INTEGER :: P2C_V(NyC) !x Lineare
1630 INTEGER :: P2C_o(NyC) !x Lineare
1631
1632 INTEGER :: P2C1_V(NyC) !x BiLineare
1633 INTEGER :: P2C2_V(NyC) !x BiLineare
1634 INTEGER :: P2C1_o(NyC) !x BiLineare
1635 INTEGER :: P2C2_o(NyC) !x BiLineare
1636
1637 REAL*8 :: diff(NrC)
1638 REAL*8 :: DEPDEP_F_WesternB(NrC)
1639 REAL*8 :: DEPDEP_F_EasternB(NrC)
1640 c----------------------------------------------------
1641 c Define CHILD model variable
1642 c----------------------------------------------------
1643 c _____________ (1) WesternB (2) EasternB
1644 c |
1645 REAL*8 :: U_F(NyC,NrC,2)
1646 REAL*8 :: V_F(NyC,NrC,2)
1647 REAL*8 :: T_F(NyC,NrC,2)
1648 REAL*8 :: S_F(NyC,NrC,2)
1649 REAL*8 :: ETA_F(NyC,NrC,2)
1650 c----------------------------------------------------
1651 PARAMETER ( terzo = 1./3.)
1652 PARAMETER (dueterzi = 2./3.)
1653 c=======================================================
1654 c (2) INTERPOLATIONS
1655 c=======================================================
1656 c (2.1) Linear for normal velocity component (u in this case)
1657 c-------------------------------------------------------
1658 DO K = 1,NrP
1659 DO J = 1,NyP-1
1660 IF (globalP4(2,J,K).eq.0.) THEN !uso la salinità come discriminante
1661 globalP1 (2,J,K) = globalP1 (2,J+1,K)
1662 ENDIF
1663 ENDDO
1664
1665 DO J = NyP,2,-1
1666 IF (globalP4(2,J,K).eq.0.) THEN !uso la salinità come discriminante
1667 globalP1 (2,J,K) = globalP1 (2,J-1,K)
1668 ENDIF
1669 ENDDO
1670 ENDDO
1671 c=======================================================
1672 c (2.1) NOT Linear but simply imposed
1673 c-------------------------------------------------------
1674 DO J = 1,3
1675 U_F(J,:,1) = globalP1(2,P2C_U(J),:)
1676 U_F(J,:,2) = globalP1(6,P2C_U(J),:)
1677 ENDDO
1678
1679 DO J = NyC-2,NyC
1680 U_F(J,:,1) = globalP1(2,P2C_U(J),:)
1681 U_F(J,:,2) = globalP1(6,P2C_U(J),:)
1682 ENDDO
1683 c=================================================
1684 DO J = 4,NyC-3,3
1685 INDC = P2C_U(J)
1686 DO K = 1,NrC
1687 c-------- WesternB ----------------------
1688 Fp = globalP1(2,INDC+1,K)
1689 Fo = globalP1(2,INDC,K)
1690 Fm = globalP1(2,INDC-1,K)
1691 c
1692 VEL_MEMO = 0.
1693 DO I = -1,1
1694 U_F(J+1+i,K,1) = ((Fp-2.*Fo+Fm)/24.)*
1695 & ((12.*float(i)**2+1.)/9.)+
1696 & ((float(i)/6.)*(Fp-Fm))+(26.*Fo-Fp-Fm)/24.
1697 VEL_MEMO = VEL_MEMO + U_F(J+1+i,K,1)
1698 ENDDO
1699
1700 VEL_MEMO = ((3.*Fo) - VEL_MEMO)/3.
1701
1702 DO I = -1,1
1703 U_F(J+1+i,K,1) = U_F(J+1+i,K,1) + VEL_MEMO
1704 ENDDO
1705
1706 c-------- EasternB ----------------------
1707 Fp = globalP1(6,INDC+1,K)
1708 Fo = globalP1(6,INDC,K)
1709 Fm = globalP1(6,INDC-1,K)
1710
1711 VEL_MEMO = 0.
1712 DO I = -1,1
1713 U_F(J+1+i,K,2) = ((Fp-2.*Fo+Fm)/24.)*
1714 & ((12.*float(i)**2+1.)/9.)+
1715 & ((float(i)/6.)*(Fp-Fm))+(26.*Fo-Fp-Fm)/24.
1716 VEL_MEMO = VEL_MEMO + U_F(J+1+i,K,2)
1717 ENDDO
1718
1719 VEL_MEMO = ((3.*Fo) - VEL_MEMO)/3.
1720
1721 DO I = -1,1
1722 U_F(J+1+i,K,2) = U_F(J+1+i,K,2) + VEL_MEMO
1723 ENDDO
1724 c------------------------------------------------------
1725 ENDDO
1726 ENDDO
1727 c-------------------------------------------------------
1728 c (2.2) BiLinear for tangent velocity component (v in this case)
1729 c-------------------------------------------------------
1730 I = 1
1731 II = WesternB
1732
1733 V_F(:,:,:) = 0.
1734
1735 DO K = 1,NrC
1736 DO J = 1,NyC
1737 x1 = Xv_P(II,P2C1_V(J))
1738 x2 = Xv_P(II,P2C2_V(J))
1739
1740 x3 = Xv_P(II+1,P2C1_V(J))
1741 x4 = Xv_P(II+1,P2C2_V(J))
1742
1743 y1 = Yv_P(II,P2C1_V(J))
1744 y2 = Yv_P(II,P2C2_V(J))
1745
1746 y3 = Yv_P(II+1,P2C1_V(J))
1747 y4 = Yv_P(II+1,P2C2_V(J))
1748
1749 x = Xv_F(I,J)
1750 y = Yv_F(I,J)
1751
1752 f1 = globalP2(1,P2C1_V(J),K)
1753 f2 = globalP2(1,P2C2_V(J),K)
1754 f3 = globalP2(2,P2C1_V(J),K)
1755 f4 = globalP2(2,P2C2_V(J),K)
1756
1757 call blint(x1,x2,x3,x4,y1,y2,y3,y4,f1,f2,f3,f4,x,y,f)
1758 V_F(J,K,1) = f
1759 ENDDO
1760 ENDDO
1761 c..............................................
1762 I = NxC
1763 II = EasternB
1764
1765 DO K = 1,NrC
1766 DO J = 1,NyC
1767 x1 = Xv_P(II,P2C1_V(J))
1768 x2 = Xv_P(II,P2C2_V(J))
1769
1770 x3 = Xv_P(II-1,P2C1_V(J))
1771 x4 = Xv_P(II-1,P2C2_V(J))
1772
1773 y1 = Yv_P(II,P2C1_V(J))
1774 y2 = Yv_P(II,P2C2_V(J))
1775
1776 y3 = Yv_P(II-1,P2C1_V(J))
1777 y4 = Yv_P(II-1,P2C2_V(J))
1778
1779 x = Xv_F(I,J)
1780 y = Yv_F(I,J)
1781
1782 f1 = globalP2(6,P2C1_V(J),K)
1783 f2 = globalP2(6,P2C2_V(J),K)
1784 f3 = globalP2(5,P2C1_V(J),K)
1785 f4 = globalP2(5,P2C2_V(J),K)
1786
1787 call blint(x1,x2,x3,x4,y1,y2,y3,y4,f1,f2,f3,f4,x,y,f)
1788 V_F(J,K,2) = f
1789 ENDDO
1790 ENDDO
1791 c-------------------------------------------------------
1792 c (2.3.1) Linear
1793 c-------------------------------------------------------
1794 c
1795 c WesternB
1796 c
1797 DO K = 1,NrP
1798 DO J = 1,NyP
1799 c
1800 DEPDEP = DEEP_P(WesternB,J,K) * DEEP_P(WesternB+1,J,K)
1801 c
1802 gammaT =(globalP3(2,J,K)-globalP3(1,J,K))*DEPDEP
1803 gammaS =(globalP4(2,J,K)-globalP4(1,J,K))*DEPDEP
1804 gammaeta =(globalP5(2,J,K)-globalP5(1,J,K))*DEPDEP
1805 cgm-------
1806 c gammaV =(globalP2(2,J,K)-globalP2(1,J,K))*DEPDEP
1807 cgm---------
1808 globalP3(1,J,K) = globalP3(1,J,K) + terzo* gammaT
1809 globalP4(1,J,K) = globalP4(1,J,K) + terzo* gammaS
1810 globalP5(1,J,K) = globalP5(1,J,K) + terzo* gammaeta
1811 cgm-------
1812 c globalP2(1,J,K) = globalP2(1,J,K) + terzo* gammaV
1813 cgm---------
1814 ENDDO
1815 c--------------------------------------------------------------
1816 DO J = 1,NyP-1
1817 IF (globalP4(1,J,K).eq.0.) THEN !uso la salinità come discriminante
1818 globalP3(1,J,K) = globalP3(1,J+1,K)
1819 globalP4(1,J,K) = globalP4(1,J+1,K)
1820 globalP5(1,J,K) = globalP5(1,J+1,K)
1821 ENDIF
1822 ENDDO
1823
1824 DO J = NyP,2,-1
1825 IF (globalP4(1,J,K).eq.0.) THEN !uso la salinità come discriminante
1826 globalP3(1,J,K) = globalP3(1,J-1,K)
1827 globalP4(1,J,K) = globalP4(1,J-1,K)
1828 globalP5(1,J,K) = globalP5(1,J-1,K)
1829 cgm---------
1830 c globalP2(1,J,K) = globalP2(1,J-1,K)
1831 cgm-------------
1832 ENDIF
1833 ENDDO
1834 c---------------------------------------------------------------
1835 ENDDO
1836
1837
1838
1839 c
1840 c EasternB
1841 c
1842 DO K = 1,NrP
1843 DO J = 1,NyP
1844 c
1845 DEPDEP = DEEP_P(EasternB,J,K) * DEEP_P(EasternB-1,J,K)
1846 c
1847 gammaT =(globalP3(5,J,K)-globalP3(6,J,K))*DEPDEP
1848 gammaS =(globalP4(5,J,K)-globalP4(6,J,K))*DEPDEP
1849 gammaeta =(globalP5(5,J,K)-globalP5(6,J,K))*DEPDEP
1850 cgm-------------
1851 c gammaV =(globalP2(5,J,K)-globalP2(6,J,K))*DEPDEP
1852 cgm----------------
1853 globalP3(6,J,K) = globalP3(6,J,K) + terzo* gammaT
1854 globalP4(6,J,K) = globalP4(6,J,K) + terzo* gammaS
1855 globalP5(6,J,K) = globalP5(6,J,K) + terzo* gammaeta
1856 cgm----------------
1857 c globalP2(6,J,K) = globalP2(6,J,K) + terzo* gammaV
1858 cgm----------------
1859 ENDDO
1860 c--------------------------------------------------------------
1861 DO J = 1,NyP-1
1862 IF (globalP4(6,J,K).eq.0.) THEN !uso la salinità come discriminante
1863 globalP3(6,J,K) = globalP3(6,J+1,K)
1864 globalP4(6,J,K) = globalP4(6,J+1,K)
1865 globalP5(6,J,K) = globalP5(6,J+1,K)
1866 cgm----------------
1867 c globalP2(6,J,K) = globalP2(6,J+1,K)
1868 cgm----------------
1869 ENDIF
1870 ENDDO
1871
1872 DO J = NyP,2,-1
1873 IF (globalP4(6,J,K).eq.0.) THEN !uso la salinità come discriminante
1874 globalP3(6,J,K) = globalP3(6,J-1,K)
1875 globalP4(6,J,K) = globalP4(6,J-1,K)
1876 globalP5(6,J,K) = globalP5(6,J-1,K)
1877 cgm----------------
1878 c globalP2(6,J,K) = globalP2(6,J-1,K)
1879 cgm----------------
1880 ENDIF
1881 ENDDO
1882 c---------------------------------------------------------------
1883 ENDDO
1884 c-------------------------------------------------------
1885 c (2.3.2) Linear
1886 c-------------------------------------------------------
1887 DO J = 1,NyC
1888 c....MASK DEEP_F
1889 IF (J.EQ.NyC) THEN !EVITO ERRORI DI CHECK BOUND x Vittorio
1890 DEPDEP_F_WesternB (:) = DEEP_F(1 ,J,:)*DEEP_F(1 ,J,:)
1891 DEPDEP_F_EasternB(:) = DEEP_F(NxC,J,:)*DEEP_F(NxC,J,:)
1892 ELSE
1893 DEPDEP_F_WesternB (:) = DEEP_F(1 ,J,:)*DEEP_F(1 ,J+1,:)
1894 DEPDEP_F_EasternB(:) = DEEP_F(NxC,J,:)*DEEP_F(NxC,J+1,:)
1895 ENDIF
1896
1897 c....WesternB...T.....................
1898 diff(:) = globalP3(1,P2C_linU(J)+1,:)-
1899 & globalP3(1,P2C_linU(J) ,:)
1900
1901 T_F(J,:,1) = globalP3(1,P2C_linU(J) ,:)+
1902 & (diff(:)/3.)*float((WO3_linU(J)))
1903 & *DEPDEP_F_WesternB(:)
1904 c.....EasternB..T.....................
1905 diff(:) = globalP3(6,P2C_linU(J)+1,:)-
1906 & globalP3(6,P2C_linU(J) ,:)
1907
1908 T_F(J,:,2) = globalP3(6,P2C_linU(J) ,:)+
1909 & (diff(:)/3.)*float((WO3_linU(J)))
1910 & *DEPDEP_F_EasternB(:)
1911 c.....WesternB...S.....................
1912 diff(:) = globalP4(1,P2C_linU(J)+1,:)-
1913 & globalP4 (1,P2C_linU(J) ,:)
1914
1915 S_F(J,:,1) = globalP4(1,P2C_linU(J) ,:)+
1916 & (diff(:)/3.)*float((WO3_linU(J)))
1917 & *DEPDEP_F_WesternB(:)
1918 c.... EasternB..S.....................
1919 diff(:) = (globalP4(6,P2C_linU(J)+1,:)-
1920 & globalP4 (6,P2C_linU(J) ,:))
1921
1922 S_F(J,:,2) = globalP4(6,P2C_linU(J) ,:)+
1923 & (diff(:)/3.)*float((WO3_linU(J)))
1924 & *DEPDEP_F_EasternB(:)
1925 c.... WesternB...Eta.....................
1926 diff(:) = globalP5(1,P2C_linU(J)+1,:)-
1927 & globalP5(1,P2C_linU(J) ,:)
1928
1929 eta_F(J,:,1) = globalP5(1,P2C_linU(J) ,:)+
1930 & (diff(:)/3.)*float((WO3_linU(J)))
1931 & *DEPDEP_F_WesternB(:)
1932 c.... EasternB..Eta.....................
1933 diff(:) = globalP5(6,P2C_linU(J)+1,:)-
1934 & globalP5(6,P2C_linU(J) ,:)
1935
1936 eta_F(J,:,2) = globalP5(6,P2C_linU(J) ,:)+
1937 & (diff(:)/3.)*float(WO3_linU(J))
1938 & *DEPDEP_F_EasternB(:)
1939
1940 cgm--------------------------
1941
1942 c.... WesternB...V.....................
1943 c diff(:) = globalP2(1,P2C_linU(J)+1,:)-
1944 c & globalP2(1,P2C_linU(J) ,:)
1945 c
1946 c V_F(J,:,1) = globalP2(1,P2C_linU(J) ,:)+
1947 c & (diff(:)/3.)*float((WO3_linU(J)))
1948 c & *DEPDEP_F_WesternB(:)
1949 c.... EasternB..V.....................
1950 c diff(:) = globalP2(6,P2C_linU(J)+1,:)-
1951 c & globalP2(6,P2C_linU(J) ,:)
1952 c
1953 c V_F(J,:,2) = globalP2(6,P2C_linU(J) ,:)+
1954 c & (diff(:)/3.)*float((WO3_linU(J)))
1955 c & *DEPDEP_F_EasternB(:)
1956 cgm----------------------------
1957 ENDDO
1958
1959
1960 8765 CONTINUE
1961 c=============== END INTERPOLAZIONI x CHILD ====================
1962
1963 RETURN
1964 END
1965
1966
1967
1968
1969
1970
1971
1972
1973
1974
1975 c================================================================
1976 SUBROUTINE BLINT(x1,x2,x3,x4,y1,y2,y3,y4,f1,f2,f3,f4,x,y,f)
1977 c================================================================
1978 C
1979 C Bilinear interpolation subroutine.
1980 C (Xi,Yi,fi) = data grid & values surounding model point (x,y)
1981 C f = interpolated value at the model grid point.
1982 IMPLICIT NONE
1983 real*4 a1,a2,a3,a4
1984 real*4 b1,b2,b3,b4
1985 real*4 x1,x2,x3,x4
1986 real*4 y1,y2,y3,y4
1987 real*4 f1,f2,f3,f4
1988 real*4 x,y,A,B,C,t,f,s
1989 C
1990 a1=x1-x2+x3-x4
1991 a2=-x1+x4
1992 a3=-x1+x2
1993 a4=x1-x
1994 b1=y1-y2+y3-y4
1995 b2=-y1+y4
1996 b3=-y1+y2
1997 b4=y1-y
1998 A=a3*b1-a1*b3
1999 B=b2*a3+b1*a4-a1*b4-a2*b3
2000 C=-a2*b4+a4*b2
2001 if(ABS(A*C).gt.0.002*B**2) then
2002 t=(-B-sqrt(B*B-4.*A*C))/(2.*A)
2003 else
2004 t=C/ABS(B)
2005 endif
2006 10 CONTINUE
2007 A=a2*b1-a1*b2
2008 B=b3*a2+b1*a4-a1*b4-a3*b2
2009 C=-a3*b4+a4*b3
2010 if(ABS(A*C).gt.0.002*B**2) then
2011 s=(-B+sqrt(B*B-4.*A*C))/(2.*A)
2012 else
2013 s=-C/ABS(B)
2014 endif
2015 20 CONTINUE
2016 f=f1*(1.-t)*(1.-s)+f2*t*(1.-s)+f3*s*t+f4*(1.-t)*s
2017 return
2018 end

  ViewVC Help
Powered by ViewVC 1.1.22