New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
modupdatebasic.F90 in vendors/AGRIF/dev/AGRIF_FILES – NEMO

source: vendors/AGRIF/dev/AGRIF_FILES/modupdatebasic.F90 @ 15504

Last change on this file since 15504 was 14975, checked in by jchanut, 3 years ago

#2638, merge new AGRIF library into trunk
  • Property svn:keywords set to Id
File size: 19.6 KB
Line 
1!
2! $Id$
3!
4!     AGRIF (Adaptive Grid Refinement In Fortran)
5!
6!     Copyright (C) 2003 Laurent Debreu (Laurent.Debreu@imag.fr)
7!                        Christophe Vouland (Christophe.Vouland@imag.fr)
8!
9!     This program is free software; you can redistribute it and/or modify
10!     it under the terms of the GNU General Public License as published by
11!     the Free Software Foundation; either version 2 of the License, or
12!     (at your option) any later version.
13!
14!     This program is distributed in the hope that it will be useful,
15!     but WITHOUT ANY WARRANTY; without even the implied warranty of
16!     MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
17!     GNU General Public License for more details.
18!
19!     You should have received a copy of the GNU General Public License
20!     along with this program; if not, write to the Free Software
21!     Foundation, Inc., 59 Temple Place -  Suite 330, Boston, MA 02111-1307, USA.
22!
23!
24!
25!> Module containing different procedures of update (copy, average, full_weighting)
26!! used in the #Agrif_Update module.
27!===================================================================================================
28!
29module Agrif_UpdateBasic
30!
31    use Agrif_Types
32
33    implicit none
34
35    integer, dimension(:,:), allocatable :: indchildcopy
36    integer, dimension(:,:), allocatable :: indchildaverage
37!
38contains
39!
40!===================================================================================================
41!  subroutine Agrif_basicupdate_copy1d
42!
43!> Carries out a copy on a parent grid (vector x) from its child grid (vector y).
44!---------------------------------------------------------------------------------------------------
45subroutine Agrif_basicupdate_copy1d ( x, y, np, nc, s_parent, s_child, ds_parent, ds_child )
46!---------------------------------------------------------------------------------------------------
47    real, dimension(np), intent(out)    :: x            !< Coarse output data to parent
48    real, dimension(nc), intent(in)     :: y            !< Fine input data from child
49    integer,             intent(in)     :: np           !< Length of parent array
50    integer,             intent(in)     :: nc           !< Length of child  array
51    real(kind=8),        intent(in)     :: s_parent     !< Parent grid position (s_root = 0)
52    real(kind=8),        intent(in)     :: s_child      !< Child  grid position (s_root = 0)
53    real(kind=8),        intent(in)     :: ds_parent    !< Parent grid dx (ds_root = 1)
54    real(kind=8),        intent(in)     :: ds_child     !< Child  grid dx (ds_root = 1)
55!---------------------------------------------------------------------------------------------------
56    integer :: i, locind_child_left, coeffraf
57!
58    coeffraf = nint(ds_parent/ds_child)
59    locind_child_left = 1 + nint((s_parent - s_child)/ds_child)
60!
61    if ( coeffraf == 1 ) then
62!CDIR ALTCODE
63        x(1:np) = y(locind_child_left:locind_child_left+np-1)
64        return
65    endif
66!
67!CDIR ALTCODE
68    do i = 1,np
69        x(i) = y(locind_child_left)
70        locind_child_left = locind_child_left + coeffraf
71    enddo
72!---------------------------------------------------------------------------------------------------
73end subroutine Agrif_basicupdate_copy1d
74!===================================================================================================
75!
76!===================================================================================================
77!  subroutine Agrif_basicupdate_copy1d_before
78!
79!> Precomputes index for a copy on a parent grid (vector x) from its child grid (vector y).
80!---------------------------------------------------------------------------------------------------
81subroutine Agrif_basicupdate_copy1d_before ( nc2, np, nc, s_parent, s_child, ds_parent, ds_child, dir )
82!---------------------------------------------------------------------------------------------------
83    integer,             intent(in)     :: nc2          !< Length of parent array
84    integer,             intent(in)     :: np           !< Length of parent array
85    integer,             intent(in)     :: nc           !< Length of child  array
86    real(kind=8),        intent(in)     :: s_parent     !< Parent grid position (s_root = 0)
87    real(kind=8),        intent(in)     :: s_child      !< Child  grid position (s_root = 0)
88    real(kind=8),        intent(in)     :: ds_parent    !< Parent grid dx (ds_root = 1)
89    real(kind=8),        intent(in)     :: ds_child     !< Child  grid dx (ds_root = 1)
90    integer,             intent(in)     :: dir          !< Direction
91!---------------------------------------------------------------------------------------------------
92    integer, dimension(:,:), allocatable    :: indchildcopy_tmp
93    integer :: i, n_old, locind_child_left, coeffraf
94!
95    coeffraf = nint(ds_parent/ds_child)
96!
97    locind_child_left = 1 + nint((s_parent - s_child)/ds_child)
98
99    if ( .not.allocated(indchildcopy) ) then
100        allocate(indchildcopy(np*nc2, 3))
101    else
102        n_old = size(indchildcopy,1)
103        if ( n_old < np*nc2 ) then
104            allocate( indchildcopy_tmp(1:n_old, 3))
105            indchildcopy_tmp = indchildcopy
106            deallocate(indchildcopy)
107            allocate(indchildcopy(np*nc2, 3))
108            indchildcopy(1:n_old,:) = indchildcopy_tmp
109            deallocate(indchildcopy_tmp)
110        endif
111    endif
112!
113    do i = 1,np
114        indchildcopy(i,dir) = locind_child_left
115        locind_child_left = locind_child_left + coeffraf
116    enddo
117!
118    do i = 2,nc2
119        indchildcopy(1+(i-1)*np:i*np,dir) = indchildcopy(1+(i-2)*np:(i-1)*np,dir) + nc
120    enddo
121!---------------------------------------------------------------------------------------------------
122end subroutine Agrif_basicupdate_copy1d_before
123!===================================================================================================
124!
125!===================================================================================================
126!  subroutine Agrif_basicupdate_copy1d_after
127!
128!> Carries out a copy on a parent grid (vector x) from its child grid (vector y)
129!! using precomputed index.
130!---------------------------------------------------------------------------------------------------
131subroutine Agrif_basicupdate_copy1d_after ( x, y, np, nc, dir )
132!---------------------------------------------------------------------------------------------------
133    real, dimension(np), intent(out)    :: x            !< Coarse output data to parent
134    real, dimension(nc), intent(in)     :: y            !< Fine input data from child
135    integer,             intent(in)     :: np           !< Length of parent array
136    integer,             intent(in)     :: nc           !< Length of child  array
137    integer,             intent(in)     :: dir          !< Direction
138!---------------------------------------------------------------------------------------------------
139    integer :: i
140!
141!CDIR ALTCODE
142    do i = 1,np
143        x(i) = y(indchildcopy(i,dir))
144    enddo
145!---------------------------------------------------------------------------------------------------
146end subroutine Agrif_basicupdate_copy1d_after
147!===================================================================================================
148!
149!===================================================================================================
150!  subroutine Agrif_basicupdate_average1d
151!
152!> Carries out an update by average on a parent grid (vector x)from its child grid (vector y).
153!---------------------------------------------------------------------------------------------------
154subroutine Agrif_basicupdate_average1d ( x, y, np, nc, s_parent, s_child, ds_parent, ds_child )
155!---------------------------------------------------------------------------------------------------
156    REAL, DIMENSION(np), intent(out)    :: x
157    REAL, DIMENSION(nc), intent(in)     :: y
158    INTEGER,             intent(in)     :: np,nc
159    REAL(kind=8),        intent(in)     :: s_parent,  s_child
160    REAL(kind=8),        intent(in)     :: ds_parent, ds_child
161!
162    INTEGER :: i, ii, locind_child_left, coeffraf
163    REAL(kind=8)    :: xpos
164    REAL ::  invcoeffraf
165    INTEGER :: nbnonnuls
166    INTEGER :: diffmod
167!
168    coeffraf = nint(ds_parent/ds_child)
169    invcoeffraf = 1./coeffraf
170!
171    if (coeffraf == 1) then
172        locind_child_left = 1 + nint((s_parent - s_child)/ds_child)
173        x(1:np) = y(locind_child_left:locind_child_left+np-1)
174        return
175    endif
176!
177    xpos = s_parent
178    x = 0.
179!
180    diffmod = 0
181!
182    IF ( mod(coeffraf,2) == 0 ) diffmod = 1
183!
184    locind_child_left = 1 + agrif_int((xpos - s_child)/ds_child)
185!
186    IF (Agrif_UseSpecialValueInUpdate) THEN
187        do i = 1,np
188            nbnonnuls = 0
189!CDIR NOVECTOR
190            do ii = -coeffraf/2+locind_child_left+diffmod, &
191                     coeffraf/2+locind_child_left
192                IF (y(ii) /= Agrif_SpecialValueFineGrid) THEN
193!                    nbnonnuls = nbnonnuls + 1
194                    x(i) = x(i) + y(ii)
195                ENDIF
196            enddo
197!            IF (nbnonnuls /= 0) THEN
198!                x(i) = x(i)/nbnonnuls
199!            ELSE
200!                x(i) = Agrif_SpecialValueFineGrid
201!            ENDIF
202            locind_child_left = locind_child_left + coeffraf
203        enddo
204    ELSE
205!
206!CDIR ALTCODE
207        do i = 1,np
208!CDIR NOVECTOR
209            do ii = -coeffraf/2+locind_child_left+diffmod, &
210                     coeffraf/2+locind_child_left
211                x(i) = x(i) + y(ii)
212            enddo
213!            x(i) = x(i)*invcoeffraf
214            locind_child_left = locind_child_left + coeffraf
215        enddo
216        IF (.not.Agrif_Update_Weights) THEN
217           x = x * invcoeffraf
218        ENDIF
219    ENDIF
220!---------------------------------------------------------------------------------------------------
221end subroutine Agrif_basicupdate_average1d
222!===================================================================================================
223
224!===================================================================================================
225!  subroutine Agrif_basicupdate_max1d
226!
227!> Carries out an update by taking the maximum on a parent grid (vector x)from its child grid (vector y).
228!---------------------------------------------------------------------------------------------------
229subroutine Agrif_basicupdate_max1d ( x, y, np, nc, s_parent, s_child, ds_parent, ds_child )
230!---------------------------------------------------------------------------------------------------
231    REAL, DIMENSION(np), intent(out)    :: x
232    REAL, DIMENSION(nc), intent(in)     :: y
233    INTEGER,             intent(in)     :: np,nc
234    REAL,                intent(in)     :: s_parent,  s_child
235    REAL,                intent(in)     :: ds_parent, ds_child
236!
237    INTEGER :: i, ii, locind_child_left, coeffraf
238    REAL    :: xpos, invcoeffraf
239    INTEGER :: nbnonnuls
240    INTEGER :: diffmod
241!
242    coeffraf = nint(ds_parent/ds_child)
243    invcoeffraf = 1./coeffraf
244!
245    if (coeffraf == 1) then
246        locind_child_left = 1 + nint((s_parent - s_child)/ds_child)
247        x(1:np) = y(locind_child_left:locind_child_left+np-1)
248        return
249    endif
250!
251    xpos = s_parent
252    x = -HUGE(1.0)
253!
254    diffmod = 0
255!
256    IF ( mod(coeffraf,2) == 0 ) diffmod = 1
257!
258    locind_child_left = 1 + agrif_int((xpos - s_child)/ds_child)
259!
260    IF (Agrif_UseSpecialValueInUpdate) THEN
261        do i = 1,np
262            nbnonnuls = 0
263!CDIR NOVECTOR
264            do ii = -coeffraf/2+locind_child_left+diffmod, &
265                     coeffraf/2+locind_child_left
266                IF (y(ii) /= Agrif_SpecialValueFineGrid) THEN
267                    x(i) = max(x(i),y(ii))
268                ENDIF
269            enddo
270            locind_child_left = locind_child_left + coeffraf
271        enddo
272    ELSE
273!
274!CDIR ALTCODE
275        do i = 1,np
276!CDIR NOVECTOR
277            do ii = -coeffraf/2+locind_child_left+diffmod, &
278                     coeffraf/2+locind_child_left
279                x(i) = max(x(i),y(ii))
280            enddo
281            locind_child_left = locind_child_left + coeffraf
282        enddo
283    ENDIF
284!---------------------------------------------------------------------------------------------------
285end subroutine Agrif_basicupdate_max1d
286!===================================================================================================
287
288!
289!===================================================================================================
290!  subroutine Average1dPrecompute
291!
292!> Carries out an update by average on a parent grid (vector x)from its child grid (vector y).
293!---------------------------------------------------------------------------------------------------
294subroutine Average1dPrecompute ( nc2, np, nc, s_parent, s_child, ds_parent, ds_child, dir )
295!---------------------------------------------------------------------------------------------------
296    INTEGER, intent(in) :: nc2, np, nc
297    REAL(kind=8),    intent(in) :: s_parent,  s_child
298    REAL(kind=8),    intent(in) :: ds_parent, ds_child
299    INTEGER, intent(in) :: dir
300!
301    INTEGER, DIMENSION(:,:), ALLOCATABLE :: indchildaverage_tmp
302    INTEGER :: i, locind_child_left, coeffraf
303    REAL(kind=8)    :: xpos
304    INTEGER :: diffmod
305!
306    coeffraf = nint(ds_parent/ds_child)
307    xpos = s_parent
308    diffmod = 0
309!
310    IF ( mod(coeffraf,2) == 0 ) diffmod = 1
311!
312    locind_child_left = 1 + agrif_int((xpos - s_child)/ds_child)
313!
314    if (.not.allocated(indchildaverage)) then
315        allocate(indchildaverage(np*nc2,3))
316    else
317        if (size(indchildaverage,1)<np*nc2) then
318            allocate( indchildaverage_tmp(size(indchildaverage,1),size(indchildaverage,2)))
319            indchildaverage_tmp = indchildaverage
320            deallocate(indchildaverage)
321            allocate(indchildaverage(np*nc2,3))
322            indchildaverage(1:size(indchildaverage_tmp,1),1:size(indchildaverage_tmp,2)) = indchildaverage_tmp
323            deallocate(indchildaverage_tmp)
324        endif
325    endif
326!
327    do i = 1,np
328        indchildaverage(i,dir)= -coeffraf/2+locind_child_left+diffmod
329        locind_child_left = locind_child_left + coeffraf
330    enddo
331!
332    do i = 2,nc2
333        indchildaverage(1+(i-1)*np:i*np,dir) = indchildaverage(1+(i-2)*np:(i-1)*np,dir) + nc
334    enddo
335!---------------------------------------------------------------------------------------------------
336end subroutine Average1dPrecompute
337!===================================================================================================
338!
339!===================================================================================================
340!  subroutine Average1dAfterCompute
341!
342!> Carries out an update by average on a parent grid (vector x) from its child grid (vector y).
343!---------------------------------------------------------------------------------------------------
344subroutine Average1dAfterCompute ( x, y, np, nc, s_parent, s_child, ds_parent, ds_child, dir )
345!---------------------------------------------------------------------------------------------------
346    REAL, DIMENSION(np), intent(inout)  :: x
347    REAL, DIMENSION(nc), intent(in)     :: y
348    INTEGER,             intent(in)     :: np, nc
349    REAL(kind=8),                intent(in)     :: s_parent,  s_child
350    REAL(kind=8),                intent(in)     :: ds_parent, ds_child
351    INTEGER,             intent(in)     :: dir
352!
353    REAL    :: invcoeffraf
354    INTEGER :: i, j, coeffraf
355    INTEGER, DIMENSION(np) :: nbnonnuls
356    REAL, DIMENSION(0:7), parameter :: invcoeff = (/1.,1.,0.5,1./3.,0.25,0.2,1./6.,1./7./)
357!
358    coeffraf = nint(ds_parent/ds_child)
359    invcoeffraf = 1./coeffraf
360!
361    IF (Agrif_UseSpecialValueInUpdate) THEN
362!
363!        nbnonnuls = 0
364        do  j = 1,coeffraf
365            do i = 1,np
366                IF (y(indchildaverage(i,dir) + j -1) /= Agrif_SpecialValueFineGrid) THEN
367!                    nbnonnuls(i) = nbnonnuls(i) + 1
368                    x(i) = x(i) +  y(indchildaverage(i,dir) + j-1 )
369                ENDIF
370            enddo
371        enddo
372        do i=1,np
373  !          x(i) = x(i)*invcoeff(nbnonnuls(i))
374  !          if (nbnonnuls(i) == 0) x(i) = Agrif_SpecialValueFineGrid
375        enddo
376!
377    ELSE
378!
379!CDIR NOLOOPCHG
380        do  j = 1,coeffraf
381!CDIR VECTOR
382            do i= 1,np
383                x(i) = x(i) + y(indchildaverage(i,dir) + j-1 )
384            enddo
385        enddo
386        IF (.not.Agrif_Update_Weights) THEN
387           x = x * invcoeffraf
388        ENDIF
389!
390    ENDIF
391
392!---------------------------------------------------------------------------------------------------
393end subroutine Average1dAfterCompute
394!===================================================================================================
395!
396!===================================================================================================
397!  subroutine Agrif_basicupdate_full_weighting1D
398!
399!> Carries out an update by full_weighting on a parent grid (vector x) from its child grid (vector y).
400!---------------------------------------------------------------------------------------------------
401subroutine Agrif_basicupdate_full_weighting1D ( x, y, np, nc, s_parent, s_child, ds_parent, ds_child )
402!---------------------------------------------------------------------------------------------------
403    real, dimension(np), intent(out)    :: x
404    real, dimension(nc), intent(in)     :: y
405    integer,             intent(in)     :: np, nc
406    real(kind=8),                intent(in)     :: s_parent,  s_child
407    real(kind=8),                intent(in)     :: ds_parent, ds_child
408!---------------------------------------------------------------------------------------------------
409    REAL(kind=8)    :: xpos, xposfin
410    INTEGER :: i, ii, diffmod
411    INTEGER :: it1, it2
412    INTEGER :: i1,  i2
413    INTEGER :: coeffraf
414    INTEGER :: locind_child_left
415    REAL    :: sumweight, invsumweight
416    REAL    :: weights(-Agrif_MaxRaff:Agrif_MaxRaff)
417
418    coeffraf = nint(ds_parent/ds_child)
419    locind_child_left = 1 + agrif_int((s_parent-s_child)/ds_child)
420!
421    if (coeffraf == 1) then
422        x(1:np) = y(locind_child_left:locind_child_left+np-1)
423        return
424    endif
425!
426    xpos = s_parent
427    x = 0.
428!
429    xposfin = s_child + ds_child * (locind_child_left - 1)
430    IF (abs(xposfin - xpos) < 0.001) THEN
431        diffmod = 0
432    ELSE
433        diffmod = 1
434    ENDIF
435!
436    if (diffmod == 1) THEN
437        invsumweight=1./(2.*coeffraf**2)
438        do i = -coeffraf,-1
439            weights(i) = invsumweight*(2*(coeffraf+i)+1)
440        enddo
441        do i = 0,coeffraf-1
442            weights(i) = weights(-(i+1))
443        enddo
444        it1 = -coeffraf
445        i1 = -(coeffraf-1)+locind_child_left
446        i2 = 2*coeffraf - 1
447       
448    else
449        invsumweight=1./coeffraf**2
450        do i = -(coeffraf-1),0
451            weights(i) = invsumweight*(coeffraf + i)
452        enddo
453        do i=1,coeffraf-1
454            weights(i) = invsumweight*(coeffraf - i)
455        enddo
456        it1 = -(coeffraf-1)
457        i1 = -(coeffraf-1)+locind_child_left
458        i2 = 2*coeffraf - 2
459
460    endif
461!
462    sumweight = 0.
463    MYLOOP : do i = 1,np
464!
465        it2 = it1
466
467!    sumweight = 0.
468   
469        do ii = i1,i1+i2
470!
471            IF (Agrif_UseSpecialValueInUpdate) THEN
472                IF (y(ii) /= Agrif_SpecialValueFineGrid) THEN
473                    x(i) = x(i) + weights(it2)*y(ii)
474!                    sumweight = sumweight+weights(it2)
475                ELSE
476                    x(i) = Agrif_SpecialValueFineGrid
477                    i1=i1+coeffraf
478                    CYCLE MYLOOP
479                ENDIF
480            ELSE
481                x(i) = x(i) + weights(it2)*y(ii)
482            ENDIF
483
484            it2 = it2+1
485!
486        enddo
487!
488        i1 = i1 + coeffraf
489!
490        enddo MYLOOP
491       
492        IF (Agrif_UseSpecialValueInUpdate) THEN
493          x = x * coeffraf ! x will be divided by coeffraf later in modupdate.F90
494        ENDIF
495       
496!---------------------------------------------------------------------------------------------------
497end subroutine Agrif_basicupdate_full_weighting1D
498!===================================================================================================
499!
500end module Agrif_UpdateBasic
Note: See TracBrowser for help on using the repository browser.