Head

GCC Code Coverage Report

Directory:	./		Exec	Total	Coverage
File:	phylmd/aaam_bud.F90	Lines:	90	92	97.8 %
Date:	2023-06-30 12:56:34	Branches:	34	38	89.5 %



! $Id: aaam_bud.F90 2350 2015-08-25 11:40:19Z emillour $

SUBROUTINE aaam_bud(iam, nlon, nlev, rjour, rsec, rea, rg, ome, plat, plon, &
    phis, dragu, liftu, phyu, dragv, liftv, phyv, p, u, v, aam, torsfc)

  USE dimphy
  USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat, klon_glo
  IMPLICIT NONE
  ! ======================================================================
  ! Auteur(s): F.Lott (LMD/CNRS) date: 20031020
  ! Object: Compute different terms of the axial AAAM Budget.
  ! No outputs, every AAM quantities are written on the IAM
  ! File.

  ! Modif : I.Musat (LMD/CNRS) date : 20041020
  ! Outputs : axial components of wind AAM "aam" and total surface torque
  ! "torsfc",
  ! but no write in the iam file.

  ! WARNING: Only valid for regular rectangular grids.
  ! REMARK: CALL DANS PHYSIQ AFTER lift_noro:
  ! CALL aaam_bud (27,klon,klev,rjourvrai,gmtime,
  ! C               ra,rg,romega,
  ! C               rlat,rlon,pphis,
  ! C               zustrdr,zustrli,zustrph,
  ! C               zvstrdr,zvstrli,zvstrph,
  ! C               paprs,u,v)

  ! ======================================================================
  ! Explicit Arguments:
  ! ==================
  ! iam-----input-I-File number where AAMs and torques are written
  ! It is a formatted file that has been opened
  ! in physiq.F
  ! nlon----input-I-Total number of horizontal points that get into physics
  ! nlev----input-I-Number of vertical levels
  ! rjour        -R-Jour compte depuis le debut de la simu (run.def)
  ! rsec         -R-Seconde de la journee
  ! rea          -R-Earth radius
  ! rg           -R-gravity constant
  ! ome          -R-Earth rotation rate
  ! plat ---input-R-Latitude en degres
  ! plon ---input-R-Longitude en degres
  ! phis ---input-R-Geopotential at the ground
  ! dragu---input-R-orodrag stress (zonal)
  ! liftu---input-R-orolift stress (zonal)
  ! phyu----input-R-Stress total de la physique (zonal)
  ! dragv---input-R-orodrag stress (Meridional)
  ! liftv---input-R-orolift stress (Meridional)
  ! phyv----input-R-Stress total de la physique (Meridional)
  ! p-------input-R-Pressure (Pa) at model half levels
  ! u-------input-R-Horizontal wind (m/s)
  ! v-------input-R-Meridional wind (m/s)
  ! aam-----output-R-Axial Wind AAM (=raam(3))
  ! torsfc--output-R-Total surface torque (=tmou(3)+tsso(3)+tbls(3))

  ! Implicit Arguments:
  ! ===================

  ! nbp_lon--common-I: Number of longitude intervals
  ! (nbp_lat-1)--common-I: Number of latitude intervals
  ! klon-common-I: Number of points seen by the physics
  ! nbp_lon*(nbp_lat-2)+2 for instance
  ! klev-common-I: Number of vertical layers
  ! ======================================================================
  ! Local Variables:
  ! ================
  ! dlat-----R: Latitude increment (Radians)
  ! dlon-----R: Longitude increment (Radians)
  ! raam  ---R: Wind AAM (3 Components, 1 & 2 Equatoriales; 3 Axiale)
  ! oaam  ---R: Mass AAM (3 Components, 1 & 2 Equatoriales; 3 Axiale)
  ! tmou-----R: Resolved Mountain torque (3 components)
  ! tsso-----R: Parameterised Moutain drag torque (3 components)
  ! tbls-----R: Parameterised Boundary layer torque (3 components)

  ! LOCAL ARRAY:
  ! ===========
  ! zs    ---R: Topographic height
  ! ps    ---R: Surface Pressure
  ! ub    ---R: Barotropic wind zonal
  ! vb    ---R: Barotropic wind meridional
  ! zlat  ---R: Latitude in radians
  ! zlon  ---R: Longitude in radians
  ! ======================================================================

  ! ARGUMENTS

  INTEGER iam, nlon, nlev
  REAL, INTENT (IN) :: rjour, rsec, rea, rg, ome
  REAL plat(nlon), plon(nlon), phis(nlon)
  REAL dragu(nlon), liftu(nlon), phyu(nlon)
  REAL dragv(nlon), liftv(nlon), phyv(nlon)
  REAL p(nlon, nlev+1), u(nlon, nlev), v(nlon, nlev)

  ! Variables locales:

  INTEGER i, j, k, l
  REAL xpi, hadley, hadday
  REAL dlat, dlon
  REAL raam(3), oaam(3), tmou(3), tsso(3), tbls(3)
  INTEGER iax
  ! IM ajout aam, torsfc
  ! aam = composante axiale du Wind AAM raam
  ! torsfc = composante axiale de (tmou+tsso+tbls)
  REAL aam, torsfc

  REAL zs(801, 401), ps(801, 401)
  REAL ub(801, 401), vb(801, 401)
  REAL ssou(801, 401), ssov(801, 401)
  REAL blsu(801, 401), blsv(801, 401)
  REAL zlon(801), zlat(401)

  CHARACTER (LEN=20) :: modname = 'aaam_bud'
  CHARACTER (LEN=80) :: abort_message



  ! PUT AAM QUANTITIES AT ZERO:

  IF (nbp_lon+1>801 .OR. nbp_lat>401) THEN
    abort_message = 'Pb de dimension dans aaam_bud'
    CALL abort_physic(modname, abort_message, 1)
  END IF

  xpi = acos(-1.)
  hadley = 1.E18
  hadday = 1.E18*24.*3600.
  IF(klon_glo.EQ.1) THEN
    dlat = xpi
  ELSE
    dlat = xpi/real(nbp_lat-1)
  ENDIF
  dlon = 2.*xpi/real(nbp_lon)

  DO iax = 1, 3
    oaam(iax) = 0.
    raam(iax) = 0.
    tmou(iax) = 0.
    tsso(iax) = 0.
    tbls(iax) = 0.
  END DO

  ! MOUNTAIN HEIGHT, PRESSURE AND BAROTROPIC WIND:

  ! North pole values (j=1):

  l = 1

  ub(1, 1) = 0.
  vb(1, 1) = 0.
  DO k = 1, nlev
    ub(1, 1) = ub(1, 1) + u(l, k)*(p(l,k)-p(l,k+1))/rg
    vb(1, 1) = vb(1, 1) + v(l, k)*(p(l,k)-p(l,k+1))/rg
  END DO

  zlat(1) = plat(l)*xpi/180.

  DO i = 1, nbp_lon + 1

    zs(i, 1) = phis(l)/rg
    ps(i, 1) = p(l, 1)
    ub(i, 1) = ub(1, 1)
    vb(i, 1) = vb(1, 1)
    ssou(i, 1) = dragu(l) + liftu(l)
    ssov(i, 1) = dragv(l) + liftv(l)
    blsu(i, 1) = phyu(l) - dragu(l) - liftu(l)
    blsv(i, 1) = phyv(l) - dragv(l) - liftv(l)

  END DO


  DO j = 2, nbp_lat-1

    ! Values at Greenwich (Periodicity)

    zs(nbp_lon+1, j) = phis(l+1)/rg
    ps(nbp_lon+1, j) = p(l+1, 1)
    ssou(nbp_lon+1, j) = dragu(l+1) + liftu(l+1)
    ssov(nbp_lon+1, j) = dragv(l+1) + liftv(l+1)
    blsu(nbp_lon+1, j) = phyu(l+1) - dragu(l+1) - liftu(l+1)
    blsv(nbp_lon+1, j) = phyv(l+1) - dragv(l+1) - liftv(l+1)
    zlon(nbp_lon+1) = -plon(l+1)*xpi/180.
    zlat(j) = plat(l+1)*xpi/180.

    ub(nbp_lon+1, j) = 0.
    vb(nbp_lon+1, j) = 0.
    DO k = 1, nlev
      ub(nbp_lon+1, j) = ub(nbp_lon+1, j) + u(l+1, k)*(p(l+1,k)-p(l+1,k+1))/rg
      vb(nbp_lon+1, j) = vb(nbp_lon+1, j) + v(l+1, k)*(p(l+1,k)-p(l+1,k+1))/rg
    END DO


    DO i = 1, nbp_lon

      l = l + 1
      zs(i, j) = phis(l)/rg
      ps(i, j) = p(l, 1)
      ssou(i, j) = dragu(l) + liftu(l)
      ssov(i, j) = dragv(l) + liftv(l)
      blsu(i, j) = phyu(l) - dragu(l) - liftu(l)
      blsv(i, j) = phyv(l) - dragv(l) - liftv(l)
      zlon(i) = plon(l)*xpi/180.

      ub(i, j) = 0.
      vb(i, j) = 0.
      DO k = 1, nlev
        ub(i, j) = ub(i, j) + u(l, k)*(p(l,k)-p(l,k+1))/rg
        vb(i, j) = vb(i, j) + v(l, k)*(p(l,k)-p(l,k+1))/rg
      END DO

    END DO

  END DO


  ! South Pole

  IF (nbp_lat-1>1) THEN
    l = l + 1
    ub(1, nbp_lat) = 0.
    vb(1, nbp_lat) = 0.
    DO k = 1, nlev
      ub(1, nbp_lat) = ub(1, nbp_lat) + u(l, k)*(p(l,k)-p(l,k+1))/rg
      vb(1, nbp_lat) = vb(1, nbp_lat) + v(l, k)*(p(l,k)-p(l,k+1))/rg
    END DO
    zlat(nbp_lat) = plat(l)*xpi/180.

    DO i = 1, nbp_lon + 1
      zs(i, nbp_lat) = phis(l)/rg
      ps(i, nbp_lat) = p(l, 1)
      ssou(i, nbp_lat) = dragu(l) + liftu(l)
      ssov(i, nbp_lat) = dragv(l) + liftv(l)
      blsu(i, nbp_lat) = phyu(l) - dragu(l) - liftu(l)
      blsv(i, nbp_lat) = phyv(l) - dragv(l) - liftv(l)
      ub(i, nbp_lat) = ub(1, nbp_lat)
      vb(i, nbp_lat) = vb(1, nbp_lat)
    END DO
  END IF


  ! MOMENT ANGULAIRE

  DO j = 1, nbp_lat-1
    DO i = 1, nbp_lon

      raam(1) = raam(1) - rea**3*dlon*dlat*0.5*(cos(zlon(i))*sin(zlat(j))*cos &
        (zlat(j))*ub(i,j)+cos(zlon(i))*sin(zlat(j+1))*cos(zlat(j+ &
        1))*ub(i,j+1)) + rea**3*dlon*dlat*0.5*(sin(zlon(i))*cos(zlat(j))*vb(i &
        ,j)+sin(zlon(i))*cos(zlat(j+1))*vb(i,j+1))

      oaam(1) = oaam(1) - ome*rea**4*dlon*dlat/rg*0.5*(cos(zlon(i))*cos(zlat( &
        j))**2*sin(zlat(j))*ps(i,j)+cos(zlon(i))*cos(zlat(j+ &
        1))**2*sin(zlat(j+1))*ps(i,j+1))

      raam(2) = raam(2) - rea**3*dlon*dlat*0.5*(sin(zlon(i))*sin(zlat(j))*cos &
        (zlat(j))*ub(i,j)+sin(zlon(i))*sin(zlat(j+1))*cos(zlat(j+ &
        1))*ub(i,j+1)) - rea**3*dlon*dlat*0.5*(cos(zlon(i))*cos(zlat(j))*vb(i &
        ,j)+cos(zlon(i))*cos(zlat(j+1))*vb(i,j+1))

      oaam(2) = oaam(2) - ome*rea**4*dlon*dlat/rg*0.5*(sin(zlon(i))*cos(zlat( &
        j))**2*sin(zlat(j))*ps(i,j)+sin(zlon(i))*cos(zlat(j+ &
        1))**2*sin(zlat(j+1))*ps(i,j+1))

      raam(3) = raam(3) + rea**3*dlon*dlat*0.5*(cos(zlat(j))**2*ub(i,j)+cos( &
        zlat(j+1))**2*ub(i,j+1))

      oaam(3) = oaam(3) + ome*rea**4*dlon*dlat/rg*0.5*(cos(zlat(j))**3*ps(i,j &
        )+cos(zlat(j+1))**3*ps(i,j+1))

    END DO
  END DO


  ! COUPLE DES MONTAGNES:


  DO j = 1, nbp_lat-1
    DO i = 1, nbp_lon
      tmou(1) = tmou(1) - rea**2*dlon*0.5*sin(zlon(i))*(zs(i,j)-zs(i,j+1))*( &
        cos(zlat(j+1))*ps(i,j+1)+cos(zlat(j))*ps(i,j))
      tmou(2) = tmou(2) + rea**2*dlon*0.5*cos(zlon(i))*(zs(i,j)-zs(i,j+1))*( &
        cos(zlat(j+1))*ps(i,j+1)+cos(zlat(j))*ps(i,j))
    END DO
  END DO

  DO j = 2, nbp_lat-1
    DO i = 1, nbp_lon
      tmou(1) = tmou(1) + rea**2*dlat*0.5*sin(zlat(j))*(zs(i+1,j)-zs(i,j))*( &
        cos(zlon(i+1))*ps(i+1,j)+cos(zlon(i))*ps(i,j))
      tmou(2) = tmou(2) + rea**2*dlat*0.5*sin(zlat(j))*(zs(i+1,j)-zs(i,j))*( &
        sin(zlon(i+1))*ps(i+1,j)+sin(zlon(i))*ps(i,j))
      tmou(3) = tmou(3) - rea**2*dlat*0.5*cos(zlat(j))*(zs(i+1,j)-zs(i,j))*( &
        ps(i+1,j)+ps(i,j))
    END DO
  END DO


  ! COUPLES DES DIFFERENTES FRICTION AU SOL:

  l = 1
  DO j = 2, nbp_lat-1
    DO i = 1, nbp_lon
      l = l + 1
      tsso(1) = tsso(1) - rea**3*cos(zlat(j))*dlon*dlat*ssou(i, j)*sin(zlat(j &
        ))*cos(zlon(i)) + rea**3*cos(zlat(j))*dlon*dlat*ssov(i, j)*sin(zlon(i &
        ))

      tsso(2) = tsso(2) - rea**3*cos(zlat(j))*dlon*dlat*ssou(i, j)*sin(zlat(j &
        ))*sin(zlon(i)) - rea**3*cos(zlat(j))*dlon*dlat*ssov(i, j)*cos(zlon(i &
        ))

      tsso(3) = tsso(3) + rea**3*cos(zlat(j))*dlon*dlat*ssou(i, j)*cos(zlat(j &
        ))

      tbls(1) = tbls(1) - rea**3*cos(zlat(j))*dlon*dlat*blsu(i, j)*sin(zlat(j &
        ))*cos(zlon(i)) + rea**3*cos(zlat(j))*dlon*dlat*blsv(i, j)*sin(zlon(i &
        ))

      tbls(2) = tbls(2) - rea**3*cos(zlat(j))*dlon*dlat*blsu(i, j)*sin(zlat(j &
        ))*sin(zlon(i)) - rea**3*cos(zlat(j))*dlon*dlat*blsv(i, j)*cos(zlon(i &
        ))

      tbls(3) = tbls(3) + rea**3*cos(zlat(j))*dlon*dlat*blsu(i, j)*cos(zlat(j &
        ))

    END DO
  END DO


  ! write(*,*) 'AAM',rsec,
  ! write(*,*) 'AAM',rjour+rsec/86400.,
  ! c      raam(3)/hadday,oaam(3)/hadday,
  ! c      tmou(3)/hadley,tsso(3)/hadley,tbls(3)/hadley

  ! write(iam,100)rjour+rsec/86400.,
  ! c      raam(1)/hadday,oaam(1)/hadday,
  ! c      tmou(1)/hadley,tsso(1)/hadley,tbls(1)/hadley,
  ! c      raam(2)/hadday,oaam(2)/hadday,
  ! c      tmou(2)/hadley,tsso(2)/hadley,tbls(2)/hadley,
  ! c      raam(3)/hadday,oaam(3)/hadday,
  ! c      tmou(3)/hadley,tsso(3)/hadley,tbls(3)/hadley
100 FORMAT (F12.5, 15(1X,F12.5))

  ! write(iam+1,*)((zs(i,j),i=1,nbp_lon),j=1,nbp_lat)
  ! write(iam+1,*)((ps(i,j),i=1,nbp_lon),j=1,nbp_lat)
  ! write(iam+1,*)((ub(i,j),i=1,nbp_lon),j=1,nbp_lat)
  ! write(iam+1,*)((vb(i,j),i=1,nbp_lon),j=1,nbp_lat)
  ! write(iam+1,*)((ssou(i,j),i=1,nbp_lon),j=1,nbp_lat)
  ! write(iam+1,*)((ssov(i,j),i=1,nbp_lon),j=1,nbp_lat)
  ! write(iam+1,*)((blsu(i,j),i=1,nbp_lon),j=1,nbp_lat)
  ! write(iam+1,*)((blsv(i,j),i=1,nbp_lon),j=1,nbp_lat)

  aam = raam(3)
  torsfc = tmou(3) + tsso(3) + tbls(3)

  RETURN
END SUBROUTINE aaam_bud


Generated by: GCOVR (Version 4.2)

Line	Branch	Exec	Source
1
2			! $Id: aaam_bud.F90 2350 2015-08-25 11:40:19Z emillour $
3
4		288	SUBROUTINE aaam_bud(iam, nlon, nlev, rjour, rsec, rea, rg, ome, plat, plon, &
5		288	phis, dragu, liftu, phyu, dragv, liftv, phyv, p, u, v, aam, torsfc)
6
7			USE dimphy
8			USE mod_grid_phy_lmdz, ONLY: nbp_lon, nbp_lat, klon_glo
9			IMPLICIT NONE
10			! ======================================================================
11			! Auteur(s): F.Lott (LMD/CNRS) date: 20031020
12			! Object: Compute different terms of the axial AAAM Budget.
13			! No outputs, every AAM quantities are written on the IAM
14			! File.
15
16			! Modif : I.Musat (LMD/CNRS) date : 20041020
17			! Outputs : axial components of wind AAM "aam" and total surface torque
18			! "torsfc",
19			! but no write in the iam file.
20
21			! WARNING: Only valid for regular rectangular grids.
22			! REMARK: CALL DANS PHYSIQ AFTER lift_noro:
23			! CALL aaam_bud (27,klon,klev,rjourvrai,gmtime,
24			! C ra,rg,romega,
25			! C rlat,rlon,pphis,
26			! C zustrdr,zustrli,zustrph,
27			! C zvstrdr,zvstrli,zvstrph,
28			! C paprs,u,v)
29
30			! ======================================================================
31			! Explicit Arguments:
32			! ==================
33			! iam-----input-I-File number where AAMs and torques are written
34			! It is a formatted file that has been opened
35			! in physiq.F
36			! nlon----input-I-Total number of horizontal points that get into physics
37			! nlev----input-I-Number of vertical levels
38			! rjour -R-Jour compte depuis le debut de la simu (run.def)
39			! rsec -R-Seconde de la journee
40			! rea -R-Earth radius
41			! rg -R-gravity constant
42			! ome -R-Earth rotation rate
43			! plat ---input-R-Latitude en degres
44			! plon ---input-R-Longitude en degres
45			! phis ---input-R-Geopotential at the ground
46			! dragu---input-R-orodrag stress (zonal)
47			! liftu---input-R-orolift stress (zonal)
48			! phyu----input-R-Stress total de la physique (zonal)
49			! dragv---input-R-orodrag stress (Meridional)
50			! liftv---input-R-orolift stress (Meridional)
51			! phyv----input-R-Stress total de la physique (Meridional)
52			! p-------input-R-Pressure (Pa) at model half levels
53			! u-------input-R-Horizontal wind (m/s)
54			! v-------input-R-Meridional wind (m/s)
55			! aam-----output-R-Axial Wind AAM (=raam(3))
56			! torsfc--output-R-Total surface torque (=tmou(3)+tsso(3)+tbls(3))
57
58			! Implicit Arguments:
59			! ===================
60
61			! nbp_lon--common-I: Number of longitude intervals
62			! (nbp_lat-1)--common-I: Number of latitude intervals
63			! klon-common-I: Number of points seen by the physics
64			! nbp_lon*(nbp_lat-2)+2 for instance
65			! klev-common-I: Number of vertical layers
66			! ======================================================================
67			! Local Variables:
68			! ================
69			! dlat-----R: Latitude increment (Radians)
70			! dlon-----R: Longitude increment (Radians)
71			! raam ---R: Wind AAM (3 Components, 1 & 2 Equatoriales; 3 Axiale)
72			! oaam ---R: Mass AAM (3 Components, 1 & 2 Equatoriales; 3 Axiale)
73			! tmou-----R: Resolved Mountain torque (3 components)
74			! tsso-----R: Parameterised Moutain drag torque (3 components)
75			! tbls-----R: Parameterised Boundary layer torque (3 components)
76
77			! LOCAL ARRAY:
78			! ===========
79			! zs ---R: Topographic height
80			! ps ---R: Surface Pressure
81			! ub ---R: Barotropic wind zonal
82			! vb ---R: Barotropic wind meridional
83			! zlat ---R: Latitude in radians
84			! zlon ---R: Longitude in radians
85			! ======================================================================
86
87			! ARGUMENTS
88
89			INTEGER iam, nlon, nlev
90			REAL, INTENT (IN) :: rjour, rsec, rea, rg, ome
91			REAL plat(nlon), plon(nlon), phis(nlon)
92			REAL dragu(nlon), liftu(nlon), phyu(nlon)
93			REAL dragv(nlon), liftv(nlon), phyv(nlon)
94			REAL p(nlon, nlev+1), u(nlon, nlev), v(nlon, nlev)
95
96			! Variables locales:
97
98			INTEGER i, j, k, l
99			REAL xpi, hadley, hadday
100			REAL dlat, dlon
101			REAL raam(3), oaam(3), tmou(3), tsso(3), tbls(3)
102			INTEGER iax
103			! IM ajout aam, torsfc
104			! aam = composante axiale du Wind AAM raam
105			! torsfc = composante axiale de (tmou+tsso+tbls)
106			REAL aam, torsfc
107
108			REAL zs(801, 401), ps(801, 401)
109			REAL ub(801, 401), vb(801, 401)
110			REAL ssou(801, 401), ssov(801, 401)
111			REAL blsu(801, 401), blsv(801, 401)
112			REAL zlon(801), zlat(401)
113
114			CHARACTER (LEN=20) :: modname = 'aaam_bud'
115			CHARACTER (LEN=80) :: abort_message
116
117
118
119			! PUT AAM QUANTITIES AT ZERO:
120
121	✓✗✗✓	288	IF (nbp_lon+1>801 .OR. nbp_lat>401) THEN
122			abort_message = 'Pb de dimension dans aaam_bud'
123			CALL abort_physic(modname, abort_message, 1)
124			END IF
125
126			xpi = acos(-1.)
127			hadley = 1.E18
128			hadday = 1.E1824.3600.
129	✓✗	288	IF(klon_glo.EQ.1) THEN
130			dlat = xpi
131			ELSE
132		288	dlat = xpi/real(nbp_lat-1)
133			ENDIF
134		288	dlon = 2.*xpi/real(nbp_lon)
135
136	✓✓	1152	DO iax = 1, 3
137		864	oaam(iax) = 0.
138		864	raam(iax) = 0.
139		864	tmou(iax) = 0.
140		864	tsso(iax) = 0.
141		1152	tbls(iax) = 0.
142			END DO
143
144			! MOUNTAIN HEIGHT, PRESSURE AND BAROTROPIC WIND:
145
146			! North pole values (j=1):
147
148			l = 1
149
150		288	ub(1, 1) = 0.
151		288	vb(1, 1) = 0.
152	✓✓	11520	DO k = 1, nlev
153		11232	ub(1, 1) = ub(1, 1) + u(l, k)*(p(l,k)-p(l,k+1))/rg
154		11520	vb(1, 1) = vb(1, 1) + v(l, k)*(p(l,k)-p(l,k+1))/rg
155			END DO
156
157		288	zlat(1) = plat(l)*xpi/180.
158
159	✓✓	9792	DO i = 1, nbp_lon + 1
160
161		9504	zs(i, 1) = phis(l)/rg
162		9504	ps(i, 1) = p(l, 1)
163		9504	ub(i, 1) = ub(1, 1)
164		9504	vb(i, 1) = vb(1, 1)
165		9504	ssou(i, 1) = dragu(l) + liftu(l)
166		9504	ssov(i, 1) = dragv(l) + liftv(l)
167		9504	blsu(i, 1) = phyu(l) - dragu(l) - liftu(l)
168		9792	blsv(i, 1) = phyv(l) - dragv(l) - liftv(l)
169
170			END DO
171
172
173	✓✓	9216	DO j = 2, nbp_lat-1
174
175			! Values at Greenwich (Periodicity)
176
177		8928	zs(nbp_lon+1, j) = phis(l+1)/rg
178		8928	ps(nbp_lon+1, j) = p(l+1, 1)
179		8928	ssou(nbp_lon+1, j) = dragu(l+1) + liftu(l+1)
180		8928	ssov(nbp_lon+1, j) = dragv(l+1) + liftv(l+1)
181		8928	blsu(nbp_lon+1, j) = phyu(l+1) - dragu(l+1) - liftu(l+1)
182		8928	blsv(nbp_lon+1, j) = phyv(l+1) - dragv(l+1) - liftv(l+1)
183			zlon(nbp_lon+1) = -plon(l+1)*xpi/180.
184		8928	zlat(j) = plat(l+1)*xpi/180.
185
186		8928	ub(nbp_lon+1, j) = 0.
187		8928	vb(nbp_lon+1, j) = 0.
188	✓✓	357120	DO k = 1, nlev
189		348192	ub(nbp_lon+1, j) = ub(nbp_lon+1, j) + u(l+1, k)*(p(l+1,k)-p(l+1,k+1))/rg
190		357120	vb(nbp_lon+1, j) = vb(nbp_lon+1, j) + v(l+1, k)*(p(l+1,k)-p(l+1,k+1))/rg
191			END DO
192
193
194	✓✓	294912	DO i = 1, nbp_lon
195
196		285696	l = l + 1
197		285696	zs(i, j) = phis(l)/rg
198		285696	ps(i, j) = p(l, 1)
199		285696	ssou(i, j) = dragu(l) + liftu(l)
200		285696	ssov(i, j) = dragv(l) + liftv(l)
201		285696	blsu(i, j) = phyu(l) - dragu(l) - liftu(l)
202		285696	blsv(i, j) = phyv(l) - dragv(l) - liftv(l)
203			zlon(i) = plon(l)*xpi/180.
204
205		285696	ub(i, j) = 0.
206		285696	vb(i, j) = 0.
207	✓✓	11436768	DO k = 1, nlev
208		11142144	ub(i, j) = ub(i, j) + u(l, k)*(p(l,k)-p(l,k+1))/rg
209		11427840	vb(i, j) = vb(i, j) + v(l, k)*(p(l,k)-p(l,k+1))/rg
210			END DO
211
212			END DO
213
214			END DO
215
216
217			! South Pole
218
219	✓✗	288	IF (nbp_lat-1>1) THEN
220		288	l = l + 1
221		288	ub(1, nbp_lat) = 0.
222		288	vb(1, nbp_lat) = 0.
223	✓✓	11520	DO k = 1, nlev
224		11232	ub(1, nbp_lat) = ub(1, nbp_lat) + u(l, k)*(p(l,k)-p(l,k+1))/rg
225		11520	vb(1, nbp_lat) = vb(1, nbp_lat) + v(l, k)*(p(l,k)-p(l,k+1))/rg
226			END DO
227		288	zlat(nbp_lat) = plat(l)*xpi/180.
228
229	✓✓	9792	DO i = 1, nbp_lon + 1
230		9504	zs(i, nbp_lat) = phis(l)/rg
231		9504	ps(i, nbp_lat) = p(l, 1)
232		9504	ssou(i, nbp_lat) = dragu(l) + liftu(l)
233		9504	ssov(i, nbp_lat) = dragv(l) + liftv(l)
234		9504	blsu(i, nbp_lat) = phyu(l) - dragu(l) - liftu(l)
235		9504	blsv(i, nbp_lat) = phyv(l) - dragv(l) - liftv(l)
236		9504	ub(i, nbp_lat) = ub(1, nbp_lat)
237		9792	vb(i, nbp_lat) = vb(1, nbp_lat)
238			END DO
239			END IF
240
241
242			! MOMENT ANGULAIRE
243
244	✓✓	9504	DO j = 1, nbp_lat-1
245	✓✓	304416	DO i = 1, nbp_lon
246
247			raam(1) = raam(1) - rea*3dlondlat0.5(cos(zlon(i))sin(zlat(j))*cos &
248			(zlat(j))ub(i,j)+cos(zlon(i))sin(zlat(j+1))*cos(zlat(j+ &
249			1))ub(i,j+1)) + rea3dlondlat0.5(sin(zlon(i))cos(zlat(j))*vb(i &
250		294912	,j)+sin(zlon(i))cos(zlat(j+1))vb(i,j+1))
251
252			oaam(1) = oaam(1) - omerea4dlondlat/rg0.5(cos(zlon(i))cos(zlat( &
253			j))*2sin(zlat(j))ps(i,j)+cos(zlon(i))cos(zlat(j+ &
254		294912	1))*2sin(zlat(j+1))*ps(i,j+1))
255
256			raam(2) = raam(2) - rea*3dlondlat0.5(sin(zlon(i))sin(zlat(j))*cos &
257			(zlat(j))ub(i,j)+sin(zlon(i))sin(zlat(j+1))*cos(zlat(j+ &
258			1))ub(i,j+1)) - rea3dlondlat0.5(cos(zlon(i))cos(zlat(j))*vb(i &
259			,j)+cos(zlon(i))cos(zlat(j+1))vb(i,j+1))
260
261			oaam(2) = oaam(2) - omerea4dlondlat/rg0.5(sin(zlon(i))cos(zlat( &
262			j))*2sin(zlat(j))ps(i,j)+sin(zlon(i))cos(zlat(j+ &
263			1))*2sin(zlat(j+1))*ps(i,j+1))
264
265			raam(3) = raam(3) + rea*3dlondlat0.5(cos(zlat(j))2ub(i,j)+cos( &
266		294912	zlat(j+1))*2ub(i,j+1))
267
268			oaam(3) = oaam(3) + omerea4dlondlat/rg0.5(cos(zlat(j))3ps(i,j &
269		9216	)+cos(zlat(j+1))*3ps(i,j+1))
270
271			END DO
272			END DO
273
274
275			! COUPLE DES MONTAGNES:
276
277
278			DO j = 1, nbp_lat-1
279		288	DO i = 1, nbp_lon
280			tmou(1) = tmou(1) - rea*2dlon0.5sin(zlon(i))(zs(i,j)-zs(i,j+1))( &
281			cos(zlat(j+1))ps(i,j+1)+cos(zlat(j))ps(i,j))
282			tmou(2) = tmou(2) + rea*2dlon0.5cos(zlon(i))(zs(i,j)-zs(i,j+1))( &
283			cos(zlat(j+1))ps(i,j+1)+cos(zlat(j))ps(i,j))
284			END DO
285			END DO
286
287	✓✓	9216	DO j = 2, nbp_lat-1
288	✓✓	294912	DO i = 1, nbp_lon
289			tmou(1) = tmou(1) + rea*2dlat0.5sin(zlat(j))(zs(i+1,j)-zs(i,j))( &
290		285696	cos(zlon(i+1))ps(i+1,j)+cos(zlon(i))ps(i,j))
291			tmou(2) = tmou(2) + rea*2dlat0.5sin(zlat(j))(zs(i+1,j)-zs(i,j))( &
292			sin(zlon(i+1))ps(i+1,j)+sin(zlon(i))ps(i,j))
293			tmou(3) = tmou(3) - rea*2dlat0.5cos(zlat(j))(zs(i+1,j)-zs(i,j))( &
294		294624	ps(i+1,j)+ps(i,j))
295			END DO
296			END DO
297
298
299			! COUPLES DES DIFFERENTES FRICTION AU SOL:
300
301			l = 1
302	✓✓	9216	DO j = 2, nbp_lat-1
303	✓✓	294912	DO i = 1, nbp_lon
304			l = l + 1
305			tsso(1) = tsso(1) - rea*3cos(zlat(j))dlondlatssou(i, j)sin(zlat(j &
306			))cos(zlon(i)) + rea3cos(zlat(j))dlondlatssov(i, j)sin(zlon(i &
307		285696	))
308
309			tsso(2) = tsso(2) - rea*3cos(zlat(j))dlondlatssou(i, j)sin(zlat(j &
310			))sin(zlon(i)) - rea3cos(zlat(j))dlondlatssov(i, j)cos(zlon(i &
311			))
312
313			tsso(3) = tsso(3) + rea*3cos(zlat(j))dlondlatssou(i, j)cos(zlat(j &
314		285696	))
315
316			tbls(1) = tbls(1) - rea*3cos(zlat(j))dlondlatblsu(i, j)sin(zlat(j &
317			))cos(zlon(i)) + rea3cos(zlat(j))dlondlatblsv(i, j)sin(zlon(i &
318		285696	))
319
320			tbls(2) = tbls(2) - rea*3cos(zlat(j))dlondlatblsu(i, j)sin(zlat(j &
321			))sin(zlon(i)) - rea3cos(zlat(j))dlondlatblsv(i, j)cos(zlon(i &
322			))
323
324			tbls(3) = tbls(3) + rea*3cos(zlat(j))dlondlatblsu(i, j)cos(zlat(j &
325		294624	))
326
327			END DO
328			END DO
329
330
331			! write(,) 'AAM',rsec,
332			! write(,) 'AAM',rjour+rsec/86400.,
333			! c raam(3)/hadday,oaam(3)/hadday,
334			! c tmou(3)/hadley,tsso(3)/hadley,tbls(3)/hadley
335
336			! write(iam,100)rjour+rsec/86400.,
337			! c raam(1)/hadday,oaam(1)/hadday,
338			! c tmou(1)/hadley,tsso(1)/hadley,tbls(1)/hadley,
339			! c raam(2)/hadday,oaam(2)/hadday,
340			! c tmou(2)/hadley,tsso(2)/hadley,tbls(2)/hadley,
341			! c raam(3)/hadday,oaam(3)/hadday,
342			! c tmou(3)/hadley,tsso(3)/hadley,tbls(3)/hadley
343			100 FORMAT (F12.5, 15(1X,F12.5))
344
345			! write(iam+1,*)((zs(i,j),i=1,nbp_lon),j=1,nbp_lat)
346			! write(iam+1,*)((ps(i,j),i=1,nbp_lon),j=1,nbp_lat)
347			! write(iam+1,*)((ub(i,j),i=1,nbp_lon),j=1,nbp_lat)
348			! write(iam+1,*)((vb(i,j),i=1,nbp_lon),j=1,nbp_lat)
349			! write(iam+1,*)((ssou(i,j),i=1,nbp_lon),j=1,nbp_lat)
350			! write(iam+1,*)((ssov(i,j),i=1,nbp_lon),j=1,nbp_lat)
351			! write(iam+1,*)((blsu(i,j),i=1,nbp_lon),j=1,nbp_lat)
352			! write(iam+1,*)((blsv(i,j),i=1,nbp_lon),j=1,nbp_lat)
353
354		288	aam = raam(3)
355		288	torsfc = tmou(3) + tsso(3) + tbls(3)
356
357		288	RETURN
358			END SUBROUTINE aaam_bud