Changeset 1459


Ignore:
Timestamp:
Jun 20, 2015, 9:22:53 AM (10 years ago)
Author:
emillour
Message:

Common dynamics: A couple of bug fixes

  • calfis[_p].F array boundaries must be explicitely specified when underlying arrays are of different sizes.
  • advect_new_p.F : missing initializations of intermediate variables at topmost layer.

EM

Location:
trunk/LMDZ.COMMON/libf
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.COMMON/libf/dyn3d/advect.F

    r1422 r1459  
    2828c   -------------
    2929
    30 #include "dimensions.h"
    31 #include "paramet.h"
    32 #include "comgeom.h"
     30      include "dimensions.h"
     31      include "paramet.h"
     32      include "comgeom.h"
    3333
    3434c   Arguments:
    3535c   ----------
    3636
    37       REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)
    38       REAL massebx(ip1jmp1,llm),masseby(ip1jm,llm),w(ip1jmp1,llm)
    39       REAL dv(ip1jm,llm),du(ip1jmp1,llm),dteta(ip1jmp1,llm)
     37      REAL,INTENT(IN) :: vcov(ip1jm,llm)
     38      REAL,INTENT(IN) :: ucov(ip1jmp1,llm)
     39      REAL,INTENT(IN) :: teta(ip1jmp1,llm)
     40      REAL,INTENT(IN) :: massebx(ip1jmp1,llm)
     41      REAL,INTENT(IN) :: masseby(ip1jm,llm)
     42      REAL,INTENT(IN) :: w(ip1jmp1,llm)
     43      REAL,INTENT(INOUT) :: dv(ip1jm,llm)
     44      REAL,INTENT(INOUT) :: du(ip1jmp1,llm)
     45      REAL,INTENT(INOUT) :: dteta(ip1jmp1,llm)
    4046
    4147c   Local:
    … …  
    5763         deuxjour = 2. * daysec
    5864
    59          DO   1  ij   = 1, ip1jmp1
    60          unsaire2(ij) = unsaire(ij) * unsaire(ij)
    61    1     CONTINUE
     65         DO ij   = 1, ip1jmp1
     66           unsaire2(ij) = unsaire(ij) * unsaire(ij)
     67         ENDDO
    6268      END IF
    6369
    … …  
    100106
    101107c
    102       DO 20 l = 1, llmm1
     108      DO l = 1, llmm1
    103109
    104110
    105111c       ......   calcul de  - w/2.    au niveau  l+1   .......
    106112
    107       DO 5  ij   = 1, ip1jmp1
    108       wsur2( ij ) = - 0.5 * w( ij,l+1 )
    109    5  CONTINUE
     113        DO ij   = 1, ip1jmp1
     114          wsur2( ij ) = - 0.5 * w( ij,l+1 )
     115        ENDDO
    110116
    111117
    112118c     .....................     calcul pour  du     ..................
    113119
    114       DO 6 ij = iip2 ,ip1jm-1
    115       ww        = wsur2 (  ij  )     + wsur2( ij+1 )
    116       uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
    117       du(ij,l)  = du(ij,l)   - ww * ( uu - uav(ij, l ) )/massebx(ij, l )
    118       du(ij,l+1)= du(ij,l+1) + ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
    119    6  CONTINUE
     120        DO ij = iip2 ,ip1jm-1
     121          ww        = wsur2 (  ij  )     + wsur2( ij+1 )
     122          uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
     123          du(ij,l)  = du(ij,l)   -ww*(uu-uav(ij,l))/massebx(ij,l)
     124          du(ij,l+1)= du(ij,l+1) +ww*(uu-uav(ij,l+1))/massebx(ij,l+1)
     125        ENDDO
    120126
    121127c     .....  correction pour  du(iip1,j,l)  ........
    … …  
    123129
    124130CDIR$ IVDEP
    125       DO   7 ij   = iip1 +iip1, ip1jm, iip1
    126       du( ij, l  ) = du( ij -iim, l  )
    127       du( ij,l+1 ) = du( ij -iim,l+1 )
    128    7  CONTINUE
     131        DO ij   = iip1 +iip1, ip1jm, iip1
     132          du( ij, l  ) = du( ij -iim, l  )
     133          du( ij,l+1 ) = du( ij -iim,l+1 )
     134        ENDDO
    129135
    130136c     .................    calcul pour   dv      .....................
    131137
    132       DO 8 ij = 1, ip1jm
    133       ww        = wsur2( ij+iip1 )   + wsur2( ij )
    134       vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
    135       dv(ij,l)  = dv(ij, l ) - ww * (vv - vav(ij, l ) )/masseby(ij, l )
    136       dv(ij,l+1)= dv(ij,l+1) + ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
    137    8  CONTINUE
     138        DO ij = 1, ip1jm
     139          ww        = wsur2( ij+iip1 )   + wsur2( ij )
     140          vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
     141          dv(ij,l)  = dv(ij, l ) - ww*(vv-vav(ij,l))/masseby(ij,l)
     142          dv(ij,l+1)= dv(ij,l+1) + ww*(vv-vav(ij,l+1))/masseby(ij,l+1)
     143        ENDDO
    138144
    139145c
    … …  
    147153c                   ...............
    148154
    149         DO 15 ij = 1, ip1jmp1
     155        DO ij = 1, ip1jmp1
    150156         ww            = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
    151157         dteta(ij, l ) = dteta(ij, l )  -  ww
    152158         dteta(ij,l+1) = dteta(ij,l+1)  +  ww
    153   15    CONTINUE
     159        ENDDO
    154160
    155       IF( conser)  THEN
    156         DO 17 ij = 1,ip1jmp1
    157         ge(ij)   = wsur2(ij) * wsur2(ij) * unsaire2(ij)
    158   17    CONTINUE
    159         gt       = SSUM( ip1jmp1,ge,1 )
    160         gtot(l)  = deuxjour * SQRT( gt/ip1jmp1 )
    161       END IF
     161        IF( conser)  THEN
     162          DO ij = 1,ip1jmp1
     163            ge(ij)   = wsur2(ij) * wsur2(ij) * unsaire2(ij)
     164          ENDDO
     165          gt       = SSUM( ip1jmp1,ge,1 )
     166          gtot(l)  = deuxjour * SQRT( gt/ip1jmp1 )
     167        END IF
    162168
    163   20  CONTINUE
     169      ENDDO ! of DO l = 1, llmm1
    164170 
    165       RETURN
    166171      END

  • trunk/LMDZ.COMMON/libf/dyn3dpar/advect_new_p.F

    r1422 r1459  
    2828c   -------------
    2929
    30 #include "dimensions.h"
    31 #include "paramet.h"
    32 #include "comgeom.h"
     30      include "dimensions.h"
     31      include "paramet.h"
     32      include "comgeom.h"
    3333
    3434c   Arguments:
    3535c   ----------
    3636
    37       REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm),teta(ip1jmp1,llm)
    38       REAL massebx(ip1jmp1,llm),masseby(ip1jm,llm),w(ip1jmp1,llm)
    39       REAL dv(ip1jm,llm),du(ip1jmp1,llm),dteta(ip1jmp1,llm)
     37      REAL,INTENT(IN) :: vcov(ip1jm,llm)
     38      REAL,INTENT(IN) :: ucov(ip1jmp1,llm)
     39      REAL,INTENT(IN) :: teta(ip1jmp1,llm)
     40      REAL,INTENT(IN) :: massebx(ip1jmp1,llm)
     41      REAL,INTENT(IN) :: masseby(ip1jm,llm)
     42      REAL,INTENT(IN) :: w(ip1jmp1,llm)
     43      REAL,INTENT(INOUT) :: dv(ip1jm,llm)
     44      REAL,INTENT(INOUT) :: du(ip1jmp1,llm)
     45      REAL,INTENT(INOUT) :: dteta(ip1jmp1,llm)
     46c   Local:
     47c   ------
     48
    4049      REAL,SAVE :: dv1(ip1jm,llm),du1(ip1jmp1,llm),dteta1(ip1jmp1,llm)
    4150      REAL,SAVE :: dv2(ip1jm,llm),du2(ip1jmp1,llm),dteta2(ip1jmp1,llm)
    42 c   Local:
    43 c   ------
    44 
    4551      REAL,SAVE :: uav(ip1jmp1,llm),vav(ip1jm,llm)
    4652      REAL wsur2(ip1jmp1)
    … …  
    6066         deuxjour = 2. * daysec
    6167
    62          DO   1  ij   = 1, ip1jmp1
     68         DO ij   = 1, ip1jmp1
    6369         unsaire2(ij) = unsaire(ij) * unsaire(ij)
    64    1     CONTINUE
     70         ENDDO
    6571      END IF
    6672
    … …  
    7783      DO ij=ijb,ije
    7884        du2(ij,1)=0.
     85        du1(ij,llm)=0.
    7986      ENDDO
    8087     
    … …  
    8592      DO ij=ijb,ije
    8693        dv2(ij,1)=0.
     94        dv1(ij,llm)=0.
    8795      ENDDO
    8896     
    … …  
    92100      DO ij=ijb,ije
    93101        dteta2(ij,1)=0.
     102        dteta1(ij,llm)=0.
    94103      ENDDO
    95104c$OMP END MASTER
    … …  
    129138           ENDDO
    130139         endif
    131          
     140
    132141      ENDDO
    133142c$OMP END DO     
    … …  
    169178
    170179c$OMP DO SCHEDULE(STATIC,OMP_CHUNK)
    171       DO 20 l = 1, llmm1
     180      DO l = 1, llmm1
    172181
    173182
    174183c       ......   calcul de  - w/2.    au niveau  l+1   .......
    175       ijb=ij_begin
    176       ije=ij_end+iip1
    177       if (pole_sud)  ije=ij_end
    178      
    179       DO 5  ij   = ijb, ije
    180       wsur2( ij ) = - 0.5 * w( ij,l+1 )
    181    5  CONTINUE
     184        ijb=ij_begin
     185        ije=ij_end+iip1
     186        if (pole_sud)  ije=ij_end
     187     
     188        DO ij   = ijb, ije
     189          wsur2( ij ) = - 0.5 * w( ij,l+1 )
     190        ENDDO
    182191
    183192
    184193c     .....................     calcul pour  du     ..................
    185194     
    186       ijb=ij_begin
    187       ije=ij_end
    188       if (pole_nord) ijb=ijb+iip1
    189       if (pole_sud)  ije=ije-iip1
    190          
    191       DO 6 ij = ijb ,ije-1
    192       ww        = wsur2 (  ij  )     + wsur2( ij+1 )
    193       uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
    194       du1(ij,l)  =  ww * ( uu - uav(ij, l ) )/massebx(ij, l )
    195       du2(ij,l+1)=  ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
    196    6  CONTINUE
     195        ijb=ij_begin
     196        ije=ij_end
     197        if (pole_nord) ijb=ijb+iip1
     198        if (pole_sud)  ije=ije-iip1
     199         
     200        DO ij = ijb ,ije-1
     201          ww        = wsur2 (  ij  )     + wsur2( ij+1 )
     202          uu        = 0.5 * ( ucov(ij,l) + ucov(ij,l+1) )
     203          du1(ij,l)  =  ww * ( uu - uav(ij, l ) )/massebx(ij, l )
     204          du2(ij,l+1)=  ww * ( uu - uav(ij,l+1) )/massebx(ij,l+1)
     205        ENDDO
    197206
    198207c     .................    calcul pour   dv      .....................
    199       ijb=ij_begin
    200       ije=ij_end
    201       if (pole_sud)  ije=ij_end-iip1
    202      
    203       DO 8 ij = ijb, ije
    204       ww        = wsur2( ij+iip1 )   + wsur2( ij )
    205       vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
    206       dv1(ij,l)  =  ww * (vv - vav(ij, l ) )/masseby(ij, l )
    207       dv2(ij,l+1)=  ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
    208    8  CONTINUE
     208        ijb=ij_begin
     209        ije=ij_end
     210        if (pole_sud)  ije=ij_end-iip1
     211     
     212        DO ij = ijb, ije
     213          ww        = wsur2( ij+iip1 )   + wsur2( ij )
     214          vv        = 0.5 * ( vcov(ij,l) + vcov(ij,l+1) )
     215          dv1(ij,l)  =  ww * (vv - vav(ij, l ) )/masseby(ij, l )
     216          dv2(ij,l+1)=  ww * (vv - vav(ij,l+1) )/masseby(ij,l+1)
     217        ENDDO
    209218
    210219c
    … …  
    220229        ije=ij_end
    221230       
    222         DO 15 ij = ijb, ije
     231        DO ij = ijb, ije
    223232         ww            = wsur2(ij) * (teta(ij,l) + teta(ij,l+1) )
    224233         dteta1(ij, l ) =   ww
    225234         dteta2(ij,l+1) =   ww
    226   15    CONTINUE
     235        ENDDO
    227236
    228237c ym ---> conser a voir plus tard
    … …  
    237246c      END IF
    238247
    239   20  CONTINUE
     248      ENDDO ! of DO l = 1, llmm1
    240249c$OMP END DO
    241250
    … …  
    279288c$OMP END DO NOWAIT     
    280289
    281       RETURN
    282290      END

  • trunk/LMDZ.COMMON/libf/dynlonlat_phylonlat/calfis.F

    r1422 r1459  
    9494c    ------------------
    9595
    96 #include "dimensions.h"
    97 #include "paramet.h"
     96      include "dimensions.h"
     97      include "paramet.h"
    9898
    9999      INTEGER ngridmx
    100100      PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm   )
    101101
    102 #include "comgeom2.h"
    103 #include "iniprint.h"
     102      include "comgeom2.h"
     103      include "iniprint.h"
    104104
    105105c    Arguments :
    … …  
    259259        call tpot2t(llm,tetamoy,tmoy,pkmoy)
    260260c SI ON TIENT COMPTE DE LA VARIATION DE G AVEC L'ALTITUDE:
    261         zlaymoy(:) = g*rad*rad/(g*rad-phimoy(:))-rad
     261        zlaymoy(1:llm) = g*rad*rad/(g*rad-phimoy(1:llm))-rad
    262262        zlevmoy(1) = phimoy(0)/g
    263263        DO l=2,llm

  • trunk/LMDZ.COMMON/libf/dynlonlat_phylonlat/calfis_p.F

    r1422 r1459  
    103103c    ------------------
    104104
    105 #include "dimensions.h"
    106 #include "paramet.h"
     105      include "dimensions.h"
     106      include "paramet.h"
    107107
    108108      INTEGER ngridmx
    109109      PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm   )
    110110
    111 #include "comgeom2.h"
    112 #include "iniprint.h"
     111      include "comgeom2.h"
     112      include "iniprint.h"
    113113#ifdef CPP_MPI
    114114      include 'mpif.h'
    … …  
    348348        call tpot2t_p(1,llm,tetamoy,tmoy,pkmoy)
    349349c SI ON TIENT COMPTE DE LA VARIATION DE G AVEC L'ALTITUDE:
    350         zlaymoy(:) = g*rad*rad/(g*rad-phimoy(:))-rad
     350        zlaymoy(1:llm) = g*rad*rad/(g*rad-phimoy(1:llm))-rad
    351351        zlevmoy(1) = phimoy(0)/g
    352352        DO l=2,llm
Note: See TracChangeset for help on using the changeset viewer.