source: LMDZ.3.3/trunk/libf/dyn3d/disvert0.F @ 1611

Last change on this file since 1611 was 2, checked in by lmdz, 25 years ago

Initial revision

  • Property svn:eol-style set to native
  • Property svn:executable set to *
  • Property svn:keywords set to Author Date Id Revision
File size: 3.4 KB
Line 
1      SUBROUTINE disvert0(pa,preff,ap,bp,dpres,presnivs,nivsigs,nivsig)
2
3c    Auteur :  P. Le Van .
4c
5      IMPLICIT NONE
6
7#include "dimensions.h"
8#include "paramet.h"
9c
10c=======================================================================
11c
12c
13c    s = sigma ** kappa   :  coordonnee  verticale
14c    dsig(l)            : epaisseur de la couche l ds la coord.  s
15c    sig(l)             : sigma a l'interface des couches l et l-1
16c    ds(l)              : distance entre les couches l et l-1 en coord.s
17c
18c=======================================================================
19c
20      REAL ap(llmp1),bp(llmp1),dpres(llm),nivsigs(llm),nivsig(llmp1)
21      REAL pa,preff,presnivs(llm)
22c
23c   declarations:
24c   -------------
25c
26      REAL sig(llm+1),dsig(llm)
27c
28      INTEGER l
29      REAL snorm
30      REAL alpha,beta,gama,delta,deltaz,h
31      INTEGER np,ierr
32      REAL pi,x
33
34      REAL SSUM
35      EXTERNAL SSUM
36c
37c-----------------------------------------------------------------------
38c
39      pi=2.*ASIN(1.)
40
41      OPEN(99,file='sigma.def',status='old',form='formatted',
42     s   iostat=ierr)
43
44c-----------------------------------------------------------------------
45c   cas 1 on lit les options dans sigma.def:
46c   ----------------------------------------
47
48      IF (ierr.eq.0) THEN
49
50      READ(99,*) deltaz
51      READ(99,*) h
52      READ(99,*) beta
53      READ(99,*) gama
54      READ(99,*) delta
55      READ(99,*) np
56      CLOSE(99)
57      alpha=deltaz/(llm*h)
58c
59
60       DO 1  l = 1, llm
61       dsig(l) = (alpha+(1.-alpha)*exp(-beta*(llm-l)))*
62     $          ( (tanh(gama*l)/tanh(gama*llm))**np +
63     $            (1.-l/FLOAT(llm))*delta )
64   1   CONTINUE
65
66       sig(1)=1.
67       DO 101 l=1,llm-1
68          sig(l+1)=sig(l)*(1.-dsig(l))/(1.+dsig(l))
69101    CONTINUE
70       sig(llm+1)=0.
71
72       DO 2  l = 1, llm
73       dsig(l) = sig(l)-sig(l+1)
74   2   CONTINUE
75c
76
77      ELSE
78c-----------------------------------------------------------------------
79c   cas 2 ancienne discretisation (LMD5...):
80c   ----------------------------------------
81
82      PRINT*,'WARNING!!! Ancienne discretisation verticale'
83
84      h=7.
85      snorm  = 0.
86      DO l = 1, llm
87         x = 2.*asin(1.) * (FLOAT(l)-0.5) / float(llm+1)
88         dsig(l) = 1.0 + 7.0 * SIN(x)**2
89         snorm = snorm + dsig(l)
90      ENDDO
91      snorm = 1./snorm
92      DO l = 1, llm
93         dsig(l) = dsig(l)*snorm
94      ENDDO
95      sig(llm+1) = 0.
96      DO l = llm, 1, -1
97         sig(l) = sig(l+1) + dsig(l)
98      ENDDO
99
100      ENDIF
101
102
103      DO l=1,llm
104        nivsigs(l) = FLOAT(l)
105      ENDDO
106
107      DO l=1,llmp1
108        nivsig(l)= FLOAT(l)
109      ENDDO
110
111c ....  Calculs  de ap(l) et de bp(l)  ....
112c .........................................
113c
114      pa        =  50 000.
115      preff     = 101 325.
116c                     .... Pascals ....
117
118      bp(llmp1) =   0.
119
120      DO l = 1, llm
121c
122ccc      ap(l) = 0.
123ccc      bp(l) = sig(l)
124
125       bp(l) = EXP( 1. -1./( sig(l)*sig(l)) )
126       ap(l) = pa * ( sig(l) - bp(l) )
127
128      ENDDO
129       ap(llmp1) = pa * ( sig(llmp1) - bp(llmp1) )
130
131      PRINT *,' BP '
132      PRINT *, bp
133      PRINT *,' AP '
134      PRINT *, ap
135
136      DO l = 1, llm
137       dpres(l) = bp(l) - bp(l+1)
138       presnivs(l) = 0.5 *( ap(l)+bp(l)*preff + ap(l+1)+bp(l+1)*preff )
139      ENDDO
140
141      PRINT *,' PRESNIVS '
142      PRINT *,presnivs
143c-----------------------------------------------------------------------
144c
145
146c-----------------------------------------------------------------------
147      RETURN
148      END
Note: See TracBrowser for help on using the repository browser.