source: LMDZ6/branches/Amaury_dev/libf/phylmd/cv3_inip.F90 @ 5128

Last change on this file since 5128 was 5116, checked in by abarral, 5 months ago

rename modules properly lmdz_*
move ismin, ismax, minmax into new lmdz_libmath.f90
(lint) uppercase fortran keywords

  • Property copyright set to
    Name of program: LMDZ
    Creation date: 1984
    Version: LMDZ5
    License: CeCILL version 2
    Holder: Laboratoire de m\'et\'eorologie dynamique, CNRS, UMR 8539
    See the license file in the root directory
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 4.4 KB
Line 
1SUBROUTINE cv3_inip()
2  ! *******************************************************************
3  ! *                                                                  *
4  ! CV3_INIP Input = choice of mixing probability laws                 *
5  !          Output = normalized coefficients of the probability laws. *
6  ! *                                                                  *
7  ! written by   : Jean-Yves Grandpeix, 06/06/2006, 19.39.27           *
8  ! modified by :                                                      *
9  ! *******************************************************************
10
11!----------------------------------------------
12!  INPUT (from Common YOMCST2 in "YOMCST2.h") :
13! iflag_mix
14! gammas   
15! alphas
16! betas
17! Fmax
18! scut
19
20!----------------------------------------------
21!  INPUT/OUTPUT (from and to Common YOMCST2 in "YOMCST2.h") :
22! qqa1
23! qqa2
24
25!----------------------------------------------
26!  OUTPUT (to Common YOMCST2 in "YOMCST2.h") :
27! Qcoef1max
28! Qcoef2max
29
30!----------------------------------------------
31
32  USE lmdz_print_control, ONLY: prt_level, lunout
33  USE lmdz_abort_physic, ONLY: abort_physic
34  IMPLICIT NONE
35
36  include "YOMCST2.h"
37
38!----------------------------------------------
39! Local variables :
40   CHARACTER (LEN=20) :: modname = 'cv3_inip'
41   CHARACTER (LEN=80) :: abort_message
42
43   REAL               :: sumcoef
44   REAL               :: sigma, aire, pdf, mu, df
45   REAL               :: ff
46
47
48  ! --   Mixing probability distribution functions
49
50  REAL qcoef1, qcoef2, qff, qfff, qmix, rmix, qmix1, rmix1, qmix2, rmix2, f
51
52  qcoef1(f) = tanh(f/gammas)
53  qcoef2(f) = (tanh(f/gammas)+gammas*log(cosh((1.-f)/gammas)/cosh(f/gammas)))
54  qff(f) = max(min(f,1.), 0.)
55  qfff(f) = min(qff(f), scut)
56  qmix1(f) = (tanh((qff(f)-fmax)/gammas)+qcoef1max)/qcoef2max
57  rmix1(f) = (gammas*log(cosh((qff(f)-fmax)/gammas))+qff(f)*qcoef1max)/ &
58    qcoef2max
59  qmix2(f) = -log(1.-qfff(f))/scut
60  rmix2(f) = (qfff(f)+(1.-qff(f))*log(1.-qfff(f)))/scut
61  qmix(f) = qqa1*qmix1(f) + qqa2*qmix2(f)
62  rmix(f) = qqa1*rmix1(f) + qqa2*rmix2(f)
63
64  ! cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
65
66
67  ! ===========================================================================
68  ! READ IN PARAMETERS FOR THE MIXING DISTRIBUTION
69  ! AND PASS THESE THROUGH A COMMON BLOCK TO SUBROUTINE CONVECT etc.
70  ! (Written by V.T.J. Phillips, 20-30/Jan/99)
71  ! ===========================================================================
72
73  ! line 1:  a flag (0 or 1) to decide whether P(F) = 1 or the general P(F)
74  ! is to be
75  ! used, followed by SCUT, which is the cut-off value of F in CONVECT
76  ! line 2:  blank
77  ! line 3:  the coefficients for the linear combination of P(F)s to
78  ! make the general P(F)
79  ! line 4:  blank
80  ! line 5:  gammas, Fmax for the cosh^2 component of P(F)
81  ! line 6:  blank
82  ! line 7:  alphas for the 1st irrational P(F)
83  ! line 8:  blank
84  ! line 9:  betas  for the 2nd irrational P(F)
85
86
87  ! c$$$        open(57,file='parameter_mix.data')
88  ! c$$$
89  ! c$$$        read(57,*) iflag_mix, scut
90  ! c$$$        read(57,*)
91  ! c$$$        IF(iflag_mix .gt. 0) THEN
92  ! c$$$              read(57,*) qqa1, qqa2
93  ! c$$$              read(57,*)
94  ! c$$$              read(57,*) gammas, Fmax
95  ! c$$$              read(57,*)
96  ! c$$$              read(57,*) alphas
97  ! c$$$         endif
98  ! c$$$         close(57)
99
100
101  IF (iflag_mix>0) THEN
102
103    ! --      Normalize Pdf weights
104
105    sumcoef = qqa1 + qqa2
106    qqa1 = qqa1/sumcoef
107    qqa2 = qqa2/sumcoef
108
109    qcoef1max = qcoef1(fmax)
110    qcoef2max = qcoef2(fmax)
111
112    sigma = 0.
113    aire = 0.0
114    pdf = 0.0
115    mu = 0.0
116    df = 0.0001
117
118    ! do ff = 0.0 + df, 1.0 - 2.*df, df
119    ff = df
120    DO WHILE (ff<=1.0-2.*df)
121      pdf = (qmix(ff+df)-qmix(ff))*(1.-ff)/df
122      aire = aire + (qmix(ff+df)-qmix(ff))*(1.-ff)
123      mu = mu + pdf*ff*df
124      IF (prt_level>9) WRITE (lunout, *) pdf, qmix(ff), aire, ff
125      ff = ff + df
126    END DO
127
128    ! do ff=0.0+df,1.0 - 2.*df,df
129    ff = df
130    DO WHILE (ff<=1.0-2.*df)
131      pdf = (qmix(ff+df)-qmix(ff))*(1.-ff)/df
132      sigma = sigma + pdf*(ff-mu)*(ff-mu)*df
133      ff = ff + df
134    END DO
135    sigma = sqrt(sigma)
136
137    IF (abs(aire-1.0)>0.02) THEN
138      WRITE (lunout, *) 'WARNING:: AREA OF MIXING PDF IS::', aire
139      abort_message = ''
140      CALL abort_physic(modname, abort_message, 1)
141    ELSE
142      PRINT *, 'Area, mean & std deviation are ::', aire, mu, sigma
143    END IF
144  END IF !  (iflag_mix .gt. 0)
145
146
147END SUBROUTINE cv3_inip
Note: See TracBrowser for help on using the repository browser.