source: LMDZ4/branches/LMDZ4_AR5/libf/phylmd/cv3_inicp.F @ 1483

Last change on this file since 1483 was 966, checked in by Laurent Fairhead, 16 years ago

Les parametres lus dans parameter_mix.data sont maintenant lus par
conf_gcm
LF

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 3.7 KB
Line 
1        SUBROUTINE cv3_inicp()
2*
3***************************************************************
4*                                                             *
5* CV3_INIP Lecture des choix de lois de probabilité de mélange*
6*          et calcul de leurs coefficients normalisés.        *
7*                                                             *
8* written by   : Jean-Yves Grandpeix, 06/06/2006, 19.39.27    *
9* modified by :                                               *
10***************************************************************
11*
12#include "YOMCST2.h"
13c
14      INTEGER iflag_clos
15c
16c --   Mixing probability distribution functions
17c
18      real Qcoef1,Qcoef2,QFF,QFFF,Qmix,Rmix,Qmix1,Rmix1,Qmix2,Rmix2,F
19      Qcoef1(F) = tanh(F/gammas)
20      Qcoef2(F) = ( tanh(F/gammas) + gammas *
21     $            log(cosh((1.- F)/gammas)/cosh(F/gammas)))
22      QFF(F) = Max(Min(F,1.),0.)
23      QFFF(F) = Min(QFF(F),scut)
24      Qmix1(F) = ( tanh((QFF(F) - Fmax)/gammas)+Qcoef1max )/
25     $           Qcoef2max
26      Rmix1(F) = ( gammas*log(cosh((QFF(F)-Fmax)/gammas))
27     1             +QFF(F)*Qcoef1max ) / Qcoef2max
28      Qmix2(F) = -Log(1.-QFFF(F))/scut
29      Rmix2(F) = (QFFF(F)+(1.-QFF(F))*Log(1.-QFFF(F)))/scut
30      Qmix(F) = qqa1*Qmix1(F) + qqa2*Qmix2(F)
31      Rmix(F) = qqa1*Rmix1(F) + qqa2*Rmix2(F)
32C
33ccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc
34c
35C
36C===========================================================================
37C       READ IN PARAMETERS FOR THE MIXING DISTRIBUTION
38C       AND PASS THESE THROUGH A COMMON BLOCK TO SUBROUTINE CONVECT etc.
39C       (Written by V.T.J. Phillips, 20-30/Jan/99)
40C===========================================================================
41C
42C   line 1:  a flag (0 or 1) to decide whether P(F) = 1 or the general P(F) is to be
43C         used, followed by SCUT, which is the cut-off value of F in CONVECT
44C   line 2:  blank
45C   line 3:  the coefficients for the linear combination of P(F)s to
46C                 make the general P(F)
47C   line 4:  blank
48C   line 5:  gammas, Fmax for the cosh^2 component of P(F)
49C   line 6:  blank
50C   line 7:  alphas for the 1st irrational P(F)
51C   line 8:  blank
52C   line 9:  betas  for the 2nd irrational P(F)
53C
54
55c        open(57,file='parameter_mix.data')
56
57c        read(57,*) iflag_clos
58c        read(57,*) iflag_mix, scut
59c        read(57,*)
60c        if(iflag_mix .gt. 0) then
61c             read(57,*) qqa1, qqa2
62c              read(57,*)
63c              read(57,*) gammas, Fmax
64c              read(57,*)
65c              read(57,*) alphas
66c         endif
67c        close(57)
68c
69      if(iflag_mix .gt. 0) then
70c
71c--      Normalize Pdf weights
72c
73        sumcoef=qqa1+qqa2
74        qqa1=qqa1/sumcoef
75        qqa2=qqa2/sumcoef
76c
77        Qcoef1max = Qcoef1(Fmax)
78        Qcoef2max = Qcoef2(Fmax)
79c
80        sigma = 0.
81        aire=0.0
82        pdf=0.0
83        mu=0.0
84        df = 0.0001
85c
86c        do ff = 0.0 + df, 1.0 - 2.*df, df
87         ff=df
88         dowhile ( ff .le. 1.0 - 2.*df )
89              pdf = (Qmix(ff+df) -  Qmix(ff)) * (1.-ff) / df
90              aire=aire+(Qmix(ff+df) - Qmix(ff)) * (1.-ff)
91              mu = mu + pdf * ff * df
92cc              write(*,*) pdf,  Qmix(ff), aire, ff
93         ff=ff+df
94         enddo
95c
96c         do ff=0.0+df,1.0 - 2.*df,df
97         ff=df
98         dowhile ( ff .le. 1.0 - 2.*df )
99              pdf = (Qmix(ff+df)- Qmix(ff)) * (1.-ff) / df
100              sigma = sigma+pdf*(ff - mu)*(ff - mu)*df
101         ff=ff+df
102         enddo
103         sigma = sqrt(sigma)
104c
105        if (abs(aire-1.0) .gt. 0.02) then
106            print *,'WARNING:: AREA OF MIXING PDF IS::', aire
107            stop
108        else
109            print *,'Area, mean & std deviation are ::', aire,mu,sigma
110        endif
111      endif     !  (iflag_mix .gt. 0)
112
113      RETURN
114      END
Note: See TracBrowser for help on using the repository browser.