The TELEMAC-MASCARET system  trunk
gaia_suspension_conv.f
Go to the documentation of this file.
1 ! ********************************
2  SUBROUTINE gaia_suspension_conv
3 ! ********************************
4  &(uconv_tel,vconv_tel,iconvf,solsys,j,flbor_w)
5 !
6 !***********************************************************************
7 ! GAIA
8 !***********************************************************************
9 !
11 ! velocities
12 !
13 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
20 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
21 !
22  USE bief
25  IMPLICIT NONE
26 !
27 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
28 !
29  TYPE(bief_obj), INTENT(IN) :: UCONV_TEL,VCONV_TEL,FLBOR_W
30  INTEGER, INTENT(IN) :: ICONVF,SOLSYS,J
31 !
32 !+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
33 !
34  DOUBLE PRECISION R1,I1,I2,A,B,AUX,LAUX,LL,USTAR,ROUSE
35  INTEGER I
36  LOGICAL GO
37  INTEGER, POINTER, DIMENSION(:) :: GLOSEG1,GLOSEG2
38 !
39 !-----------------------------------------------------------------------
40 !
41 ! OPTVF : TENS 0 : NORMAL
42 ! 1 : ADVECTION FIELD DOES NOT SATISFY
43 ! CONTINUITY
44 !
45 ! OPTVF : UNITS 0 : CONSTANT = 0
46 ! 1 : CHI-TUAN CONSTANT
47 ! 2 : LEO POSTMA CONSTANT
48 ! SEE CVTRVF IN BIEF AND
49 ! V5.7 RELEASE NOTES
50 !
51  solsys_gai=solsys
54  gloseg1=>mesh%GLOSEG%I(1:mesh%GLOSEG%DIM1)
55  gloseg2=>mesh%GLOSEG%I(mesh%GLOSEG%DIM1+1:2*mesh%GLOSEG%DIM1)
56 !
57  yaflulim_gai=.false.
58  go=.false.
59 !
60  IF(corr_conv.AND.(.NOT.sedco(j))) THEN
61 !
62  ll=log(30.d0)
63 !
64  DO i = 1, npoin
65 !
66  IF(tob%R(i).GT.zero) THEN
67 !
68 ! Here we use TOBCW_MEAN instead of TOB to define USTAR.
69 ! In case there are no wave effects in the simulation,
70 ! TOBCW_MEAN = TOB. In case there are waves, the
71 ! superimposed effects of waves and current are taken
72 ! into account to establish the suspension concentration
73 ! profile through USTAR.
74 
75  ustar = sqrt(tobcw_mean%R(i)/xmve)
76 !
77 ! B --> KS/H
78 !
79 ! AUX = 1.D0 + KARMAN*SQRT(2.D0/MAX(CF%R(I),ZERO))
80 ! B = 30.D0*EXP(-AUX)
81 !
82  b = ksr%R(i) /max(hn%R(i),1.1d0*ksr%R(i))
83  a = zref%R(i)/max(hn%R(i),1.1d0*zref%R(i))
84 !
85 ! TAKES MAX VALUE OF A = ZREF/H AND B=KSR/H
86  a=max(a,b)
87 !
88 ! SIMPLIFIED VERSION
89  rouse=min(xwc(j)/max(ustar,zero),1.d0)/karman
90  r1= 1.d0-rouse
91  laux=log(a)
92 !
93  IF(abs(r1).LT.1.d-8) THEN
94  i1= -laux
95  i2= -laux**2/2.d0
96  ELSE
97  aux=a**r1
98  i1=(1.d0-aux)/r1
99  i2=-(i1+laux*aux)/r1
100  ENDIF
101 !
102 ! AUX=LOG(A/30.D0)
103  aux=laux - ll
104  t1%R(i)=-(i2-aux*i1)/(i1*(aux+1.d0))
105 !
106  ELSE
107 !
108  t1%R(i)=1.d0
109 !
110  ENDIF
111 !
112 ! CHECKS 0
113 !
114  t1%R(i)=min(t1%R(i),1.d0)
115  t1%R(i)=max(t1%R(i),0.d0)
116 !
117  ENDDO
118 !
119 ! DEPENDING ON ADVECTION SCHEME : LIMITATION OF VELOCITY OR FLUXES
120 !
121  IF(iconvf.EQ.13.OR.iconvf.EQ.14) THEN
122 !
123 ! LIMITATION OF FLUXES WITH FLULIM
124 !
125  uconv_gai%R=>uconv_tel%R
126  vconv_gai%R=>vconv_tel%R
127  DO i=1,mesh%NSEG
128  flulim_gai%R(i)=0.5d0*(t1%R(gloseg1(i))+t1%R(gloseg2(i)))
129  ENDDO
130  yaflulim_gai=.true.
131 !
132  ELSE
133 !
134 ! LIMITATION OF VELOCITY
135 !
136  solsys_gai=1
137  DO i = 1,npoin
138  uconv_gai%R(i) = t1%R(i)*u2d%R(i)
139  vconv_gai%R(i) = t1%R(i)*v2d%R(i)
140  ENDDO
141  yaflulim_gai=.false.
142 !
143  ENDIF
144 !
145 ! ADVECTION FORM WHICH ACCEPTS AN ADVECTION FIELD
146 ! THAT DOES NOT SATISFY CONTINUITY + LEO-POSTMA CONSTANT
147 !
148 ! WITH 12: MASS CONSERVATION BUT NO MONOTONICITY
149 ! THE CORRECT THEORY
150  optvf_gai=12
151 !
152 ! WITH 2: MONOTONICITY BUT NO MASS CONSERVATION
153 ! WRONG THEORY
154 ! OPTVF=2
155 !
156 ! OPTVF=2 IS POSSIBLE BUT WITH MASS CONSERVATION SPOILED
157 ! THE UNIT (HERE 2) IS REDONE IN CVDFTR ACCORDING TO THE
158 ! VALUE OF RESOL, SO IT IS NOT IMPORTANT HERE.
159 !
160  CALL os('X=Y ',x=flbor_gai,y=flbor_w)
161  CALL osbd('X=CXY ',flbor_gai,t1,t1,1.d0,mesh)
162 !
163  ELSE
164 !
165 ! POINTERS ARE USED TO AVOID COPY
166 !
167 ! HERE UCONV_TEL IS PASSED ON
168  uconv_gai%R=>uconv_tel%R
169  vconv_gai%R=>vconv_tel%R
170 ! ADVECTION FORM THAT REQUIRES AN ADVECTION FIELD
171 ! THAT SATISFIES CONTINUITY + LEO-POSTMA CONSTANT
172  optvf_gai=2
173 !
174  ENDIF
175 !
176  RETURN
177  END
type(bief_obj), target u2d
Components of depth-averaged velocity.
type(bief_obj), target uconv_gai
Components of velocity vectors.
type(bief_obj), target tob
Bed shear stress [n/m2].
integer solsys_gai
Choose the advection field in cvdftr.
logical corr_conv
Correction on convection velocity.
type(bief_obj), target v2d
double precision, dimension(:), pointer x
2d coordinates of the mesh
double precision, dimension(nsiclm), target xwc
Settling velocities.
double precision xmve
Water density (from steering file of T2D or T3D)
double precision, dimension(:), pointer save_uconv
Save the velocity fields in the suspension computation.
type(bief_obj), target flbor_gai
Flux at the boundaries.
double precision, dimension(:), pointer y
type(bief_obj), target flulim_gai
Flux limitation per segment.
double precision zero
Parameter used for clipping variables or testing values against zero.
double precision, dimension(:), pointer save_vconv
subroutine osbd(OP, X, Y, Z, C, MESH)
Definition: osbd.f:7
type(bief_obj), target hn
Water depth.
logical, dimension(nsiclm) sedco
Cohesive sediments (for each class)
type(bief_obj), target zref
Reference elevation.
type(bief_obj), pointer t1
Aliases for work vectors in tb.
type(bief_obj), target vconv_gai
subroutine os(OP, X, Y, Z, C, IOPT, INFINI, ZERO)
Definition: os.f:7
type(bief_mesh), target mesh
Mesh structure.
type(bief_obj), target tobcw_mean
Mean of total current + wave shear stress.
subroutine gaia_suspension_conv(UCONV_TEL, VCONV_TEL, ICONVF, SOLSYS, J, FLBOR_W)
double precision karman
Karman constant.
type(bief_obj), target ksr
Ripple bed roughness.
integer optvf_gai
Option for finite volumes (see cvtrvf)
logical yaflulim_gai
Logical for modification of boundary fluxes.
integer, pointer npoin
Number of 2d points in the mesh.
Definition: bief.f:3