Ignore:
Timestamp:
Oct 21, 2024, 2:58:45 PM (23 hours ago)
Author:
abarral
Message:

Convert fixed-form to free-form sources .F -> .{f,F}90
(WIP: some .F remain, will be handled in subsequent commits)

File:
1 moved

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/dyn3dmem/advect_new_loc.F90

    r5245 r5246  
    22! $Header$
    33!
    4       SUBROUTINE advect_new_loc(ucov,vcov,teta,w,massebx,masseby,
    5      &                        du,dv,dteta)
    6       USE parallel_lmdz
    7       USE write_field_loc
    8       USE advect_new_mod
    9       USE comconst_mod, ONLY: daysec
    10       USE logic_mod, ONLY: conser
    11      
    12       IMPLICIT NONE
    13 c=======================================================================
    14 c
    15 c   Auteurs:  P. Le Van , Fr. Hourdin  .
    16 c   -------
    17 c
    18 c   Objet:
    19 c   ------
    20 c
    21 c   *************************************************************
    22 c   .... calcul des termes d'advection vertic.pour u,v,teta,q ...
    23 c   *************************************************************
    24 c        ces termes sont ajoutes a du,dv,dteta et dq .
    25 c  Modif F.Forget 03/94 : on retire q de advect
    26 c
    27 c=======================================================================
    28 c-----------------------------------------------------------------------
    29 c   Declarations:
    30 c   -------------
    31 
    32       include "dimensions.h"
    33       include "paramet.h"
    34       include "comgeom.h"
    35 
    36 c   Arguments:
    37 c   ----------
    38 
    39       REAL vcov(ijb_v:ije_v,llm),ucov(ijb_u:ije_u,llm)
    40       REAL teta(ijb_u:ije_u,llm)
    41       REAL massebx(ijb_u:ije_u,llm),masseby(ijb_v:ije_v,llm)
    42       REAL w(ijb_u:ije_u,llm)
    43       REAL dv(ijb_v:ije_v,llm),du(ijb_u:ije_u,llm)
    44       REAL dteta(ijb_u:ije_u,llm)
    45 c   Local:
    46 c   ------
    47 
    48       REAL wsur2(ijb_u:ije_u)
    49       REAL unsaire2(ijb_u:ije_u), ge(ijb_u:ije_u)
    50       REAL deuxjour, ww, gt, uu, vv
    51 
    52       INTEGER ij,l,ijb,ije
    53       EXTERNAL  SSUM
    54       REAL      SSUM
    55 
    56 
    57 
    58 c-----------------------------------------------------------------------
    59 c   2. Calculs preliminaires:
    60 c   -------------------------
    61      
    62       IF (conser.AND.1==0)  THEN
    63          deuxjour = 2. * daysec
    64 
    65          DO   1  ij   = 1, ip1jmp1
    66          unsaire2(ij) = unsaire(ij) * unsaire(ij)
    67    1     CONTINUE
    68       END IF
    69 
    70 
    71 c------------------  -yy ----------------------------------------------
    72 c   .  Calcul de     u
    73 
    74 c$OMP MASTER
    75       ijb=ij_begin
    76       ije=ij_end
    77       if (pole_nord) ijb=ijb+iip1
    78       if (pole_sud)  ije=ije-iip1
    79 
    80       DO ij=ijb,ije
    81         du2(ij,1)=0.
    82         du1(ij,llm)=0.
    83       ENDDO
    84      
    85       ijb=ij_begin
    86       ije=ij_end
    87       if (pole_sud)  ije=ij_end-iip1
    88      
    89       DO ij=ijb,ije
    90         dv2(ij,1)=0.
    91         dv1(ij,llm)=0.
    92       ENDDO
    93      
    94       ijb=ij_begin
    95       ije=ij_end
    96 
    97       DO ij=ijb,ije
    98         dteta2(ij,1)=0.
    99         dteta1(ij,llm)=0.
    100       ENDDO
    101 c$OMP END MASTER
    102 c$OMP BARRIER
    103  
    104 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
    105       DO  l=1,llm
    106          
    107          ijb=ij_begin
    108          ije=ij_end
    109          if (pole_nord) ijb=ijb+iip1
    110          if (pole_sud)  ije=ije-iip1
    111          
    112 c        DO    ij     = iip2, ip1jmp1
    113 c            uav(ij,l) = 0.25 * ( ucov(ij,l) + ucov(ij-iip1,l) )
    114 c        ENDDO
    115 
    116 c        DO    ij     = iip2, ip1jm
    117 c            uav(ij,l) = uav(ij,l) + uav(ij+iip1,l)
    118 c        ENDDO
    119          
    120          DO    ij     = ijb, ije
    121                  
    122            uav(ij,l)=0.25*(ucov(ij,l)+ucov(ij-iip1,l))
    123      .               +0.25*(ucov(ij+iip1,l)+ucov(ij,l))
    124          ENDDO
    125          
    126          if (pole_nord) then
    127            DO      ij         = 1, iip1
    128               uav(ij      ,l) = 0.
    129            ENDDO
    130          endif
    131          
    132          if (pole_sud) then
    133            DO      ij         = 1, iip1
    134               uav(ip1jm+ij,l) = 0.
    135            ENDDO
    136          endif
    137          
    138       ENDDO
    139 c$OMP END DO     
    140 c      call write_field3d_p('uav',reshape(uav,(/iip1,jjp1,llm/)))
    141      
    142 c------------------  -xx ----------------------------------------------
    143 c   .  Calcul de     v
    144      
    145       ijb=ij_begin
    146       ije=ij_end
    147       if (pole_sud)  ije=ij_end-iip1
    148 
    149 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
    150       DO  l=1,llm
    151          
    152          DO    ij   = ijb+1, ije
    153            vav(ij,l) = 0.25 * ( vcov(ij,l) + vcov(ij-1,l) )
    154          ENDDO
    155          
    156          DO    ij   = ijb,ije,iip1
    157           vav(ij,l) = vav(ij+iim,l)
    158          ENDDO
    159          
    160          
    161          DO    ij   = ijb, ije-1
    162           vav(ij,l) = vav(ij,l) + vav(ij+1,l)
    163          ENDDO
    164          
    165          DO    ij       = ijb, ije, iip1
    166           vav(ij+iim,l) = vav(ij,l)
    167          ENDDO
    168          
    169       ENDDO
    170 c$OMP END DO
    171 c      call write_field3d_p('vav',reshape(vav,(/iip1,jjm,llm/)))
    172 
    173 c-----------------------------------------------------------------------
    174 c$OMP BARRIER
    175 
    176 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    177       DO 20 l = 1, llmm1
    178 
    179 
    180 c      ......   calcul de  - w/2.    au niveau  l+1   .......
    181       ijb=ij_begin
    182       ije=ij_end+iip1
    183       if (pole_sud)  ije=ij_end
    184      
    185       DO 5   ij   = ijb, ije
    186       wsur2( ij ) = - 0.5 * w( ij,l+1 )
    187    5  CONTINUE
    188 
    189 
    190 c    .....................     calcul pour  du     ..................
    191      
    192       ijb=ij_begin
    193       ije=ij_end
    194       if (pole_nord) ijb=ijb+iip1
    195       if (pole_sud)  ije=ije-iip1
    196          
    197       DO 6 ij = ijb ,ije-1
    198       ww        = wsur2 (  ij  )     + wsur2( ij+1 )
    199       uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
    200       du1(ij,l)  =  ww * ( uu - uav(ij, l ) )/massebx(ij, l )
    201       du2(ij,l+1)=  ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
    202    6  CONTINUE
    203 
    204 c    .................    calcul pour   dv      .....................
    205       ijb=ij_begin
    206       ije=ij_end
    207       if (pole_sud)  ije=ij_end-iip1
    208      
    209       DO 8 ij = ijb, ije
    210       ww        = wsur2( ij+iip1 )   + wsur2( ij )
    211       vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
    212       dv1(ij,l)  =  ww * (vv - vav(ij, l ) )/masseby(ij, l )
    213       dv2(ij,l+1)=  ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
    214    8  CONTINUE
    215 
    216 c
    217 
    218 c    ............................................................
    219 c    ...............    calcul pour   dh      ...................
    220 c    ............................................................
    221 
    222 c                       ---z
    223 c       calcul de  - d( teta  * w )      qu'on ajoute a   dh
    224 c                   ...............
    225         ijb=ij_begin
    226         ije=ij_end
    227        
    228         DO 15 ij = ijb, ije
    229          ww            = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
    230          dteta1(ij, l ) =   ww
    231          dteta2(ij,l+1) =   ww
    232   15    CONTINUE
    233 
    234 c ym ---> conser a voir plus tard
    235 
    236 c      IF( conser)  THEN
    237 c       
    238 c        DO 17 ij = 1,ip1jmp1
    239 c        ge(ij)   = wsur2(ij) * wsur2(ij) * unsaire2(ij)
    240 c  17    CONTINUE
    241 c        gt       = SSUM( ip1jmp1,ge,1 )
    242 c        gtot(l)  = deuxjour * SQRT( gt/ip1jmp1 )
    243 c      END IF
    244 
    245   20  CONTINUE
    246 c$OMP END DO
    247 
    248       ijb=ij_begin
    249       ije=ij_end
    250       if (pole_nord) ijb=ijb+iip1
    251       if (pole_sud)  ije=ije-iip1
    252 #ifdef DEBUG_IO   
    253        CALL WriteField_u('du_bis',du)     
     4SUBROUTINE advect_new_loc(ucov,vcov,teta,w,massebx,masseby, &
     5        du,dv,dteta)
     6  USE parallel_lmdz
     7  USE write_field_loc
     8  USE advect_new_mod
     9  USE comconst_mod, ONLY: daysec
     10  USE logic_mod, ONLY: conser
     11
     12  IMPLICIT NONE
     13  !=======================================================================
     14  !
     15  !   Auteurs:  P. Le Van , Fr. Hourdin  .
     16  !   -------
     17  !
     18  !   Objet:
     19  !   ------
     20  !
     21  !   *************************************************************
     22  !   .... calcul des termes d'advection vertic.pour u,v,teta,q ...
     23  !   *************************************************************
     24  !    ces termes sont ajoutes a du,dv,dteta et dq .
     25  !  Modif F.Forget 03/94 : on retire q de advect
     26  !
     27  !=======================================================================
     28  !-----------------------------------------------------------------------
     29  !   Declarations:
     30  !   -------------
     31
     32  include "dimensions.h"
     33  include "paramet.h"
     34  include "comgeom.h"
     35
     36  !   Arguments:
     37  !   ----------
     38
     39  REAL :: vcov(ijb_v:ije_v,llm),ucov(ijb_u:ije_u,llm)
     40  REAL :: teta(ijb_u:ije_u,llm)
     41  REAL :: massebx(ijb_u:ije_u,llm),masseby(ijb_v:ije_v,llm)
     42  REAL :: w(ijb_u:ije_u,llm)
     43  REAL :: dv(ijb_v:ije_v,llm),du(ijb_u:ije_u,llm)
     44  REAL :: dteta(ijb_u:ije_u,llm)
     45  !   Local:
     46  !   ------
     47
     48  REAL :: wsur2(ijb_u:ije_u)
     49  REAL :: unsaire2(ijb_u:ije_u), ge(ijb_u:ije_u)
     50  REAL :: deuxjour, ww, gt, uu, vv
     51
     52  INTEGER :: ij,l,ijb,ije
     53  EXTERNAL  SSUM
     54  REAL :: SSUM
     55
     56
     57
     58  !-----------------------------------------------------------------------
     59  !   2. Calculs preliminaires:
     60  !   -------------------------
     61
     62  IF (conser.AND.1==0)  THEN
     63     deuxjour = 2. * daysec
     64
     65     DO  ij   = 1, ip1jmp1
     66     unsaire2(ij) = unsaire(ij) * unsaire(ij)
     67     END DO
     68  END IF
     69
     70
     71  !------------------  -yy ----------------------------------------------
     72  !   .  Calcul de     u
     73
     74!$OMP MASTER
     75  ijb=ij_begin
     76  ije=ij_end
     77  if (pole_nord) ijb=ijb+iip1
     78  if (pole_sud)  ije=ije-iip1
     79
     80  DO ij=ijb,ije
     81    du2(ij,1)=0.
     82    du1(ij,llm)=0.
     83  ENDDO
     84
     85  ijb=ij_begin
     86  ije=ij_end
     87  if (pole_sud)  ije=ij_end-iip1
     88
     89  DO ij=ijb,ije
     90    dv2(ij,1)=0.
     91    dv1(ij,llm)=0.
     92  ENDDO
     93
     94  ijb=ij_begin
     95  ije=ij_end
     96
     97  DO ij=ijb,ije
     98    dteta2(ij,1)=0.
     99    dteta1(ij,llm)=0.
     100  ENDDO
     101!$OMP END MASTER
     102!$OMP BARRIER
     103
     104!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     105  DO  l=1,llm
     106
     107     ijb=ij_begin
     108     ije=ij_end
     109     if (pole_nord) ijb=ijb+iip1
     110     if (pole_sud)  ije=ije-iip1
     111
     112      ! DO    ij     = iip2, ip1jmp1
     113      !    uav(ij,l) = 0.25 * ( ucov(ij,l) + ucov(ij-iip1,l) )
     114      ! ENDDO
     115
     116      ! DO    ij     = iip2, ip1jm
     117      !    uav(ij,l) = uav(ij,l) + uav(ij+iip1,l)
     118      ! ENDDO
     119
     120     DO    ij     = ijb, ije
     121
     122       uav(ij,l)=0.25*(ucov(ij,l)+ucov(ij-iip1,l)) &
     123             +0.25*(ucov(ij+iip1,l)+ucov(ij,l))
     124     ENDDO
     125
     126     if (pole_nord) then
     127       DO      ij         = 1, iip1
     128          uav(ij      ,l) = 0.
     129       ENDDO
     130     endif
     131
     132     if (pole_sud) then
     133       DO      ij         = 1, iip1
     134          uav(ip1jm+ij,l) = 0.
     135       ENDDO
     136     endif
     137
     138  ENDDO
     139!$OMP END DO
     140   ! call write_field3d_p('uav',reshape(uav,(/iip1,jjp1,llm/)))
     141
     142  !------------------  -xx ----------------------------------------------
     143  !   .  Calcul de     v
     144
     145  ijb=ij_begin
     146  ije=ij_end
     147  if (pole_sud)  ije=ij_end-iip1
     148
     149!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     150  DO  l=1,llm
     151
     152     DO    ij   = ijb+1, ije
     153       vav(ij,l) = 0.25 * ( vcov(ij,l) + vcov(ij-1,l) )
     154     ENDDO
     155
     156     DO    ij   = ijb,ije,iip1
     157      vav(ij,l) = vav(ij+iim,l)
     158     ENDDO
     159
     160
     161     DO    ij   = ijb, ije-1
     162      vav(ij,l) = vav(ij,l) + vav(ij+1,l)
     163     ENDDO
     164
     165     DO    ij       = ijb, ije, iip1
     166      vav(ij+iim,l) = vav(ij,l)
     167     ENDDO
     168
     169  ENDDO
     170!$OMP END DO
     171    ! call write_field3d_p('vav',reshape(vav,(/iip1,jjm,llm/)))
     172
     173  !-----------------------------------------------------------------------
     174!$OMP BARRIER
     175
     176!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     177  DO l = 1, llmm1
     178
     179
     180    ! ......   calcul de  - w/2.    au niveau  l+1   .......
     181  ijb=ij_begin
     182  ije=ij_end+iip1
     183  if (pole_sud)  ije=ij_end
     184
     185  DO   ij   = ijb, ije
     186  wsur2( ij ) = - 0.5 * w( ij,l+1 )
     187  END DO
     188
     189
     190  ! .....................     calcul pour  du     ..................
     191
     192  ijb=ij_begin
     193  ije=ij_end
     194  if (pole_nord) ijb=ijb+iip1
     195  if (pole_sud)  ije=ije-iip1
     196
     197  DO ij = ijb ,ije-1
     198  ww        = wsur2 (  ij  )     + wsur2( ij+1 )
     199  uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
     200  du1(ij,l)  =  ww * ( uu - uav(ij, l ) )/massebx(ij, l )
     201  du2(ij,l+1)=  ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
     202  END DO
     203
     204  ! .................    calcul pour   dv      .....................
     205  ijb=ij_begin
     206  ije=ij_end
     207  if (pole_sud)  ije=ij_end-iip1
     208
     209  DO ij = ijb, ije
     210  ww        = wsur2( ij+iip1 )   + wsur2( ij )
     211  vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
     212  dv1(ij,l)  =  ww * (vv - vav(ij, l ) )/masseby(ij, l )
     213  dv2(ij,l+1)=  ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
     214  END DO
     215
     216  !
     217
     218  ! ............................................................
     219  ! ...............    calcul pour   dh      ...................
     220  ! ............................................................
     221
     222  !                   ---z
     223  !   calcul de  - d( teta  * w )      qu'on ajoute a   dh
     224  !               ...............
     225    ijb=ij_begin
     226    ije=ij_end
     227
     228    DO ij = ijb, ije
     229     ww            = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
     230     dteta1(ij, l ) =   ww
     231     dteta2(ij,l+1) =   ww
     232    END DO
     233
     234  ! ym ---> conser a voir plus tard
     235
     236   ! IF( conser)  THEN
     237  !
     238  !    DO 17 ij = 1,ip1jmp1
     239  !    ge(ij)   = wsur2(ij) * wsur2(ij) * unsaire2(ij)
     240  !  17    CONTINUE
     241  !    gt       = SSUM( ip1jmp1,ge,1 )
     242  !    gtot(l)  = deuxjour * SQRT( gt/ip1jmp1 )
     243  !  END IF
     244
     245  END DO
     246!$OMP END DO
     247
     248  ijb=ij_begin
     249  ije=ij_end
     250  if (pole_nord) ijb=ijb+iip1
     251  if (pole_sud)  ije=ije-iip1
     252#ifdef DEBUG_IO
     253   CALL WriteField_u('du_bis',du)
    254254#endif
    255 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)       
    256       DO l=1,llm
    257         DO ij=ijb,ije-1
    258           du(ij,l)=du(ij,l)+du2(ij,l)-du1(ij,l)
    259         ENDDO
    260 
    261         DO   ij   = ijb+iip1-1, ije, iip1
    262          du( ij, l  ) = du( ij -iim, l  )
    263         ENDDO
    264       ENDDO
    265 c$OMP END DO NOWAIT
    266 #ifdef DEBUG_IO   
    267       CALL WriteField_u('du1',du1)     
    268       CALL WriteField_u('du2',du2)     
    269       CALL WriteField_u('du_bis',du)     
     255!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     256  DO l=1,llm
     257    DO ij=ijb,ije-1
     258      du(ij,l)=du(ij,l)+du2(ij,l)-du1(ij,l)
     259    ENDDO
     260
     261    DO   ij   = ijb+iip1-1, ije, iip1
     262     du( ij, l  ) = du( ij -iim, l  )
     263    ENDDO
     264  ENDDO
     265!$OMP END DO NOWAIT
     266#ifdef DEBUG_IO
     267  CALL WriteField_u('du1',du1)
     268  CALL WriteField_u('du2',du2)
     269  CALL WriteField_u('du_bis',du)
    270270#endif
    271       ijb=ij_begin
    272       ije=ij_end
    273       if (pole_sud)  ije=ij_end-iip1
    274 
    275 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)     
    276       DO l=1,llm
    277         DO ij=ijb,ije
    278           dv(ij,l)=dv(ij,l)+dv2(ij,l)-dv1(ij,l)
    279         ENDDO
    280       ENDDO
    281 c$OMP END DO NOWAIT     
    282       ijb=ij_begin
    283       ije=ij_end
    284 
    285 c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)   
    286       DO l=1,llm
    287         DO ij=ijb,ije
    288           dteta(ij,l)=dteta(ij,l)+dteta2(ij,l)-dteta1(ij,l)
    289         ENDDO
    290       ENDDO
    291 c$OMP END DO NOWAIT     
    292 
    293       RETURN
    294       END
     271  ijb=ij_begin
     272  ije=ij_end
     273  if (pole_sud)  ije=ij_end-iip1
     274
     275!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     276  DO l=1,llm
     277    DO ij=ijb,ije
     278      dv(ij,l)=dv(ij,l)+dv2(ij,l)-dv1(ij,l)
     279    ENDDO
     280  ENDDO
     281!$OMP END DO NOWAIT
     282  ijb=ij_begin
     283  ije=ij_end
     284
     285!$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
     286  DO l=1,llm
     287    DO ij=ijb,ije
     288      dteta(ij,l)=dteta(ij,l)+dteta2(ij,l)-dteta1(ij,l)
     289    ENDDO
     290  ENDDO
     291!$OMP END DO NOWAIT
     292
     293  RETURN
     294END SUBROUTINE advect_new_loc
Note: See TracChangeset for help on using the changeset viewer.