The TELEMAC-MASCARET system  trunk
maxslope_gaia.f
Go to the documentation of this file.
1 ! ************************
2  SUBROUTINE maxslope_gaia
3 ! ************************
4 !
5 !***********************************************************************
6 ! GAIA
7 !***********************************************************************
8 !
11 !
12 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
14 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
15 !
16  USE interface_gaia, ex_maxslope => maxslope_gaia
17  USE bief
19  USE interface_parallel, ONLY: p_max
21  IMPLICIT NONE
22 !
23 !!-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
24 !
25 
26 !
27 !!-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
28 !
29  INTEGER IELEM,I1,I2,I3,I,IG1,IG2,IR1,IR2
30  INTEGER IPOIN,ILAYER,ISAND,IMUD,ICLA
31  DOUBLE PRECISION X2,X3,Y2,Y3,Z2,Z3,A,B,L,ZC,DEUXSURF,TANSL
32  DOUBLE PRECISION QGLI(3),QG1,QG2,QR1,QR2
33  DOUBLE PRECISION FAC1,FAC2,FAC3
34  DOUBLE PRECISION QELEM(nelem,3)
35  DOUBLE PRECISION POURCENTAGE_EPAI_TRANSFERT
36  DOUBLE PRECISION EPAI_GLISS,EPAI_GLISS_LAYER,FLUX_LOC
37  DOUBLE PRECISION TEST_GLISS
38  INTEGER NB_PDT_GLISS ! passer en mot clé ?
39  TYPE(bief_obj),POINTER :: EVOL
40 !
41  LOGICAL CASE2
42 !
43  INTRINSIC sqrt,min,max,tan
44 !
45 !-----------------------------------------------------------------------
46 !
47  evol => tb%ADR(1)%P ! WORK ARRAY
48 !
49  pi = 4.d0 * atan( 1.d0 )
50  tansl = tan( pi*phised/180.d0 )
52  nb_pdt_gliss = 10
53 !
54  CALL cpstvc(unsv2d,evol)
55 !
56 !======================================================================
57 ! TOTAL MASS OF SEDIMENT IN EACH POINT BEFORE SLIDING
58 ! JUST FOR CHECK THAT EVOLUTION OF MASS = 0
60 ! INITIALISES
61  CALL os('X=0 ',x=evol_mm)
62 ! EVOL_MM%R() CONTAINED -MASS OF SEDIMENT BEFORE SLIDING
63 ! LATER WE ADD MASS OF SEDIMENT AFTER SLIDING TO HAVE EVOLUTION
64 !
65  IF(nsand.GE.1) THEN
66  DO ipoin = 1,npoin
67  DO ilayer = 1,nomblay
68  DO isand = 1,nsand
69  evol_mm%R(ipoin) = evol_mm%R(ipoin)
70  & - mass_sand(isand,ilayer,ipoin)
71  ENDDO
72  ENDDO
73  ENDDO
74  ENDIF
75  IF(nmud.GE.1) THEN
76  DO ipoin = 1,npoin
77  DO ilayer = 1,nomblay
78  DO imud = 1,nmud
79  evol_mm%R(ipoin) = evol_mm%R(ipoin)
80  & - mass_mud(imud,ilayer,ipoin)
81  ENDDO
82  ENDDO
83  ENDDO
84  ENDIF
85 !
86 !-------------------------------------------------------------------------
87 !
88  DO ipoin = 1,npoin
89  IF(nmud.GE.1) THEN
90  DO ilayer = 1,nomblay
91  DO imud = 1,nmud
92  flux_mass_mud(imud,ilayer,ipoin)=0.d0
93  ENDDO
94  ENDDO
95  ENDIF
96  IF(nsand.GE.1) THEN
97  DO ilayer = 1,nomblay
98  DO isand = 1,nsand
99  flux_mass_sand(isand,ilayer,ipoin)=0.d0
100  ENDDO
101  ENDDO
102  ENDIF
103  ENDDO
104 !
105 ! SLIDING IS BASED ON ZF BUT IS APPLIED LAYER BY LAYER.
106 ! SLIDING IS FIRST APPLIED ON THE FIRST LAYER, A NEW ZF IS COMPUTED.
107 ! THEN, IF NECESSARY, SLIDING IS APPLIED ON THE OTHER LAYERS.
108  DO ilayer=1,nomblay
109 !
110  CALL os('X=0 ',x=evol)
111 !
112 ! TO CHECK IF ANY SLIDING IS DONE AND SHORTCUT IF NOT
113  test_gliss = 0.d0
114 !
115 ! FIRST LOOP ON ELEMENTS. THIS FIRST LOOP ENABLES TO COMPUTE THE
116 ! TRANSFERS OF VOLUMES BETWEEN THE 3 POINTS INSIDE EACH ELEMENT,
117 ! AND TO COMPUTE THE MAXIMUM POSSIBLE EROSION FOR EACH POINTS,
118 ! IN ORDER TO ENABLE LIMITATION OF FLUXES ON RIGID BEDS.
119 ! THE ACTUAL LIMITATION OF FLUXES ON RIGID BEDS, AND THE APPLYING
120 ! OF FLUXES TO MASSES OF THE LAYERS ARES DONE IN THE SECOND LOOP
121 ! ON ELEMENTS.
122  DO ielem=1,nelem
123 !
124  i1=mesh%IKLE%I(ielem)
125  i2=mesh%IKLE%I(ielem+nelem)
126  i3=mesh%IKLE%I(ielem+2*nelem)
127 !
128  x2=mesh%XEL%R(ielem+nelem)
129  x3=mesh%XEL%R(ielem+2*nelem)
130  y2=mesh%YEL%R(ielem+nelem)
131  y3=mesh%YEL%R(ielem+2*nelem)
132  z2=zf%R(i2)-zf%R(i1)
133  z3=zf%R(i3)-zf%R(i1)
134 !
135 ! TWICE THE TRIANGLE AREA
136 !
137  deuxsurf=x2*y3-x3*y2
138 !
139 ! AVERAGE BOTTOM IN THE ELEMENT
140 !
141  zc=(zf%R(i1)+zf%R(i2)+zf%R(i3))/3.d0
142 !
143 ! COMPONENTS OF BOTTOM GRADIENT
144 !
145  a=(z2*y3-z3*y2)/deuxsurf
146  b=(z3*x2-z2*x3)/deuxsurf
147 !
148 ! CORRECTING FACTOR ON SLOPE
149 !
150  l=min(1.d0,tansl/max(sqrt(a**2+b**2),1.d-8))
151 !
152 ! HERE THE EVOLUTIONS ARE MULTIPLIED BY SURFAC/3
153 ! BECAUSE THE REAL EVOLUTION TAKING INTO ACCOUNT OTHER ELEMENTS
154 ! WILL NEED A FACTOR (SURFAC/3)/(INTEGRAL OF BASIS)
155 !
156  qelem(ielem,1)=(1.d0/nb_pdt_gliss)*
157  & (1.d0-l)*(zc-zf%R(i1))*deuxsurf/6.d0
158  qelem(ielem,2)=(1.d0/nb_pdt_gliss)*
159  & (1.d0-l)*(zc-zf%R(i2))*deuxsurf/6.d0
160  qelem(ielem,3)=(1.d0/nb_pdt_gliss)*
161  & (1.d0-l)*(zc-zf%R(i3))*deuxsurf/6.d0
162 !!
163 ! QELEM AND EVOL ARE VOLUMES.
164 ! EVOL IS USED ONLY FOR FLUX LIMITATION ON RIGID BEDS
165 ! ONLY QELEM CORRESPONDING TO EROSIONS ARE ADDED TO EVOL (PESSIMISTIC),
166 ! OTHERWISE THE LIMITATION IS NOT ENOUGH
167 !
168  IF(qelem(ielem,1).LE.0.d0) THEN
169  evol%R(i1)=evol%R(i1)+qelem(ielem,1)
170  ENDIF
171  IF(qelem(ielem,2).LE.0.d0) THEN
172  evol%R(i2)=evol%R(i2)+qelem(ielem,2)
173  ENDIF
174  IF(qelem(ielem,3).LE.0.d0) THEN
175  evol%R(i3)=evol%R(i3)+qelem(ielem,3)
176  ENDIF
177 !
178  ENDDO ! END OF FIRST LOOP ON ELEMENTS
179 !
180 ! ADDING EVOL IN PARALLEL
181  IF(ncsize.GT.1) THEN
182  CALL parcom(evol,2,mesh)
183  ENDIF
184 !
185  DO i=1,npoin
186  test_gliss=test_gliss+abs(evol%R(i))
187  ENDDO
188 !
189 ! GET PARALLEL MAX TO MAKE SURE EACH PARTITION DOES THE SAME (BECAUSE THER ARE MORE PARCOM CALLS BELOW)
190  test_gliss = p_max(test_gliss)
191 !
192 ! IF NO SLIDING NEEDED, END OF SUBROUTINE
193  IF(test_gliss.LE.1.d-8) GOTO 100
194 !
195 
196 !MDS refaire une boucle sur les éléments pour appliquer le facteur ES(ILAYER)/EVOL aux evolutions
197 !MDS et pour appliquer les Q en terme de masse.
198 !
199 ! SECOND LOOP ON ELEMENTS
200 ! HERE THE LIMITATION OF FLUXES ON RIGID BED IS MADE, AND
201 ! THE FLUXES (QELEM, VOLUMETRIC) ARE CONVERTED IN FLUXES OF MASS
202  DO ielem=1,nelem
203 !
204  i1=mesh%IKLE%I(ielem)
205  i2=mesh%IKLE%I(ielem+nelem)
206  i3=mesh%IKLE%I(ielem+2*nelem)
207 !
208 ! LIMITATION FACTOR LESS THAN 1 WHEN LIMITATION OVER RIGID BED
209 
210  fac1=1.d0
211  IF(qelem(ielem,1).LE.0.d0) THEN
212  fac1=max(0.d0,
213  & min(1.d0,
214  & es(i1,ilayer)/max(-evol%R(i1)*unsv2d%R(i1),1.d-15)))
215  ENDIF
216 !
217  fac2=1.d0
218  IF(qelem(ielem,2).LE.0.d0) THEN
219  fac2=max(0.d0,
220  & min(1.d0,
221  & es(i2,ilayer)/max(-evol%R(i2)*unsv2d%R(i2),1.d-15)))
222  ENDIF
223 !
224  fac3=1.d0
225  IF(qelem(ielem,3).LE.0.d0) THEN
226  fac3=max(0.d0,
227  & min(1.d0,
228  & es(i3,ilayer)/max(-evol%R(i3)*unsv2d%R(i3),1.d-15)))
229  ENDIF
230 !
231  qgli(1)=qelem(ielem,1)*fac1*fac2*fac3
232  qgli(2)=qelem(ielem,2)*fac1*fac2*fac3
233  qgli(3)=qelem(ielem,3)*fac1*fac2*fac3
234 !
235 ! IG1 AND IG2 : POINTS THAT GIVE
236 ! IR1 AND IR2 : POINTS THAT RECEIVE
237 ! CASE2: TWO POINTS GIVE TO THE THIRD ONE (THE OTHER CASE IS
238 ! ONE POINT GIVES TO THE TWO OTHERS)
239  case2=.false.
240 !
241 ! PARAMETERISING TO REDUCE THE 6 CASES TO 2
242 !
243 !! Q ARE VOLUMES
244 ! FLUX LEAVING IS NEGATIVE.
245  IF(qgli(1).GT.0.d0) THEN
246  IF(qgli(2).GT.0.d0) THEN
247 ! 3 GIVES TO 1 AND 2
248  ig1=i3
249  qg1=qgli(3)
250  ir1=i1
251  qr1=qgli(1)
252  ir2=i2
253  qr2=qgli(2)
254  ELSE
255  IF(qgli(3).GT.0.d0) THEN
256 ! 2 GIVES TO 1 AND 3
257  ig1=i2
258  qg1=qgli(2)
259  ir1=i1
260  qr1=qgli(1)
261  ir2=i3
262  qr2=qgli(3)
263  ELSE
264 ! 2 AND 3 GIVE TO 1
265  ig1=i2
266  qg1=qgli(2)
267  ig2=i3
268  qg2=qgli(3)
269  ir1=i1
270  qr1=qgli(1)
271  case2=.true.
272  ENDIF
273  ENDIF
274  ELSE
275  IF(qgli(2).GT.0.d0) THEN
276  IF(qgli(3).GT.0.d0) THEN
277 ! 1 GIVES TO 2 AND 3
278  ig1=i1
279  qg1=qgli(1)
280  ir1=i2
281  qr1=qgli(2)
282  ir2=i3
283  qr2=qgli(3)
284  ELSE
285 ! 1 AND 3 GIVE TO 2
286  ig1=i1
287  qg1=qgli(1)
288  ig2=i3
289  qg2=qgli(3)
290  ir1=i2
291  qr1=qgli(2)
292  case2=.true.
293  ENDIF
294  ELSE
295 ! 1 AND 2 GIVE TO 3
296  ig1=i1
297  qg1=qgli(1)
298  ig2=i2
299  qg2=qgli(2)
300  ir1=i3
301  qr1=qgli(3)
302  case2=.true.
303  ENDIF
304  ENDIF
305 !
306 !
307  IF(case2) THEN
308 !
309 ! THE TWO DONORS CASE : IG1 AND IG2 GIVE TO IR1
310 ! FILLS FLUX_MASS_SAND AND FLUX_MASS_MUD
311 ! WHICH ARE IN DECLARATION_GAIA BUT USED ONLY LOCALLY.
312 ! WE DO NOT NEED ANY INFO ON QR1.
313 
314 ! FIRST TRANSFER OF MASS : FROM IG1 TO IR1
315  epai_gliss = -qg1*unsv2d%R(ig1)
316  IF(epai_gliss.GT.0.d0) THEN
317 ! IF EPAI_GLISS > ES THE WHOLE LAYER WILL SLIDE. OTHERWISE ONLY A PERCENTAGE OF THE LAYER
318  pourcentage_epai_transfert=
319  & max(0.d0,min(1.d0,epai_gliss/es(ig1,ilayer)))
320 !
321  IF(nmud.GE.1) THEN
322  DO imud = 1,nmud
323  flux_loc=-min(max(0.d0,pourcentage_epai_transfert*
324  & mass_mud(imud,ilayer,ig1)),
325  & mass_mud(imud,ilayer,ig1))
326  flux_mass_mud(imud,ilayer,ig1)=
327  & flux_mass_mud(imud,ilayer,ig1)
328  & +flux_loc
329  flux_mass_mud(imud,ilayer,ir1)=
330  & flux_mass_mud(imud,ilayer,ir1)
331  & -flux_loc*v2dpar%R(ig1)/v2dpar%R(ir1)
332  ENDDO
333  ENDIF
334  IF(nsand.GE.1) THEN
335  DO isand = 1,nsand
336  flux_loc=-min(max(0.d0,pourcentage_epai_transfert*
337  & mass_sand(isand,ilayer,ig1)),
338  & mass_sand(isand,ilayer,ig1))
339  flux_mass_sand(isand,ilayer,ig1)=
340  & flux_mass_sand(isand,ilayer,ig1)
341  & +flux_loc
342  flux_mass_sand(isand,ilayer,ir1)=
343  & flux_mass_sand(isand,ilayer,ir1)
344  & -flux_loc*v2dpar%R(ig1)/v2dpar%R(ir1)
345  ENDDO
346  ENDIF
347  ENDIF ! ENDIF EPAI_GLISS > 0
348 ! END OF TRANSFER FROM IG1 TO IR1
349 !
350 !!! 2) FROM IG2 TO IR1
351  epai_gliss = -qg2*unsv2d%R(ig2)
352  IF(epai_gliss.GT.0.d0) THEN
353  pourcentage_epai_transfert =
354  & max(0.d0,min(1.d0,epai_gliss/es(ig2,ilayer)))
355  epai_gliss_layer =
356  & pourcentage_epai_transfert*es(ig2,ilayer)
357 !
358  IF(nmud.GE.1) THEN
359  DO imud = 1,nmud
360  flux_loc=-min(max(0.d0,pourcentage_epai_transfert*
361  & mass_mud(imud,ilayer,ig2)),
362  & mass_mud(imud,ilayer,ig2))
363  flux_mass_mud(imud,ilayer,ig2) =
364  & flux_mass_mud(imud,ilayer,ig2)
365  & +flux_loc
366  flux_mass_mud(imud,ilayer,ir1) =
367  & flux_mass_mud(imud,ilayer,ir1)
368  & -flux_loc*v2dpar%R(ig2)/v2dpar%R(ir1)
369  ENDDO
370  ENDIF
371  IF(nsand.GE.1) THEN
372  DO isand = 1,nsand
373  flux_loc=-min(max(0.d0,pourcentage_epai_transfert*
374  & mass_sand(isand,ilayer,ig2)),
375  & mass_sand(isand,ilayer,ig2))
376  flux_mass_sand(isand,ilayer,ig2)=
377  & flux_mass_sand(isand,ilayer,ig2)
378  & +flux_loc
379  flux_mass_sand(isand,ilayer,ir1)=
380  & flux_mass_sand(isand,ilayer,ir1)
381  & -flux_loc*v2dpar%R(ig2)/v2dpar%R(ir1)
382  ENDDO
383  ENDIF
384  ENDIF ! ENDIF EPAI_GLISS > 0
385 ! END OF TRANSFER FROM IG1 TO IR2
386 !
387  ELSE
388 !
389 ! THE ONE DONOR CASE : IG1 GIVES TO IR1 AND IR2
390 ! QR1 AND QR2 ONLY USED TO COMPUTE THE RATIO TO KNOW
391 ! THE REPARTITION OF QG1 BETWWEN IR1 AND IR2
392  epai_gliss = -qg1*unsv2d%R(ig1)
393  IF(epai_gliss.GT.0.d0) THEN
394  pourcentage_epai_transfert =
395  & max(0.d0,min(1.d0,epai_gliss/es(ig1,ilayer)))
396  epai_gliss_layer =
397  & pourcentage_epai_transfert*es(ig1,ilayer)
398  IF(nmud.GE.1) THEN
399  DO imud = 1,nmud
400  flux_loc=-min(max(0.d0,pourcentage_epai_transfert*
401  & mass_mud(imud,ilayer,ig1)),
402  & mass_mud(imud,ilayer,ig1))
403  flux_mass_mud(imud,ilayer,ig1) =
404  & flux_mass_mud(imud,ilayer,ig1)
405  & +flux_loc
406  flux_mass_mud(imud,ilayer,ir1) =
407  & flux_mass_mud(imud,ilayer,ir1)
408  & -qr1/(qr1+qr2)*flux_loc*v2dpar%R(ig1)/v2dpar%R(ir1)
409  flux_mass_mud(imud,ilayer,ir2) =
410  & flux_mass_mud(imud,ilayer,ir2)
411  & -qr2/(qr1+qr2)*flux_loc*v2dpar%R(ig1)/v2dpar%R(ir2)
412  ENDDO
413  ENDIF
414  IF(nsand.GE.1) THEN
415  DO isand = 1,nsand
416  flux_loc=-min(max(0.d0,pourcentage_epai_transfert*
417  & mass_sand(isand,ilayer,ig1)),
418  & mass_sand(isand,ilayer,ig1))
419  flux_mass_sand(isand,ilayer,ig1) =
420  & flux_mass_sand(isand,ilayer,ig1)
421  & +flux_loc
422  flux_mass_sand(isand,ilayer,ir1) =
423  & flux_mass_sand(isand,ilayer,ir1)
424  & -qr1/(qr1+qr2)*flux_loc*v2dpar%R(ig1)/v2dpar%R(ir1)
425  flux_mass_sand(isand,ilayer,ir2) =
426  & flux_mass_sand(isand,ilayer,ir2)
427  & -qr2/(qr1+qr2)*flux_loc*v2dpar%R(ig1)/v2dpar%R(ir2)
428  ENDDO
429  ENDIF
430  ENDIF ! ENDIF EPAI_GLISS > 0
431 ! END OF CASE IG1 > IR1 AND IR2
432 !
433  ENDIF ! CASE2
434 !
435  ENDDO ! END OF SECOND LOOP ON ELEMENTS
436 !
437 ! ADDING VOLUMES IN PARALLEL
438  IF(ncsize.GT.1) THEN
439  IF(nmud.GE.1) THEN
440  DO imud=1,nmud
441  DO ipoin = 1,npoin
442  t2%R(ipoin) = flux_mass_mud(imud,ilayer,ipoin)
443  ENDDO
444  CALL parcom(t2,2,mesh)
445  DO ipoin = 1,npoin
446  flux_mass_mud(imud,ilayer,ipoin) = t2%R(ipoin)
447  ENDDO
448  ENDDO
449  ENDIF
450  IF(nsand.GE.1) THEN
451  DO isand=1,nsand
452  DO ipoin = 1,npoin
453  t2%R(ipoin) = flux_mass_sand(isand,ilayer,ipoin)
454  ENDDO
455  CALL parcom(t2,2,mesh)
456  DO ipoin = 1,npoin
457  flux_mass_sand(isand,ilayer,ipoin) = t2%R(ipoin)
458  ENDDO
459  ENDDO
460  ENDIF
461  ENDIF
462 !
463 ! UPDATING THE MASSES
464 !
465  DO ipoin = 1,npoin
466  IF(nmud.GE.1) THEN
467  DO imud = 1,nmud
468  mass_mud(imud,ilayer,ipoin) =
469  & mass_mud(imud,ilayer,ipoin)
470  & +flux_mass_mud(imud,ilayer,ipoin)
471  ENDDO
472  ENDIF
473  IF(nsand.GE.1) THEN
474  DO isand = 1,nsand
475  mass_sand(isand,ilayer,ipoin) =
476  & mass_sand(isand,ilayer,ipoin)
477  & +flux_mass_sand(isand,ilayer,ipoin)
478  ENDDO
479  ENDIF
480  ENDDO
481 !
482  IF(vsmtype==0) THEN
483  IF(debug.GT.0) WRITE(lu,*)'BED1_UPDATE'
484  CALL bed1_update(zr,zf,volu2d)
485  IF(debug.GT.0) WRITE(lu,*)'END BED1_UPDATE'
486  ELSEIF(vsmtype==1) THEN
487 ! CVSM
488 ! !==== conversion of mass to thickness and fraction =============!
489 ! EVCL_MB and MASS_SAND given in [ kg/m**2 ]
490  DO icla = 1,nsicla ! convert for layer-1 for each sand class,
491  DO i = 1,npoin
492  zfcl_c%ADR(icla)%P%R(i)
493  & = evcl_mb%ADR(icla)%P%R(i) * mpa2t(icla) ! thickness evolution
494  IF(es(i,1).GT.0.d0) THEN
495  avail(i,1,icla)
496  & = mass_sand(icla,1,i) * mpa2t(icla) / es(i,1) ! mass to fraction
497  ELSE
498  avail(i,1,icla) = 0.d0
499  ENDIF
500  zf%R(i) = zf%R(i)+zfcl_c%ADR(icla)%P%R(i) ! new ZF
501  ENDDO
502  ENDDO
503  IF(debug.GT.0) WRITE(lu,*)'CVSP_MAIN_GAIA'
505  IF(debug.GT.0) WRITE(lu,*)'CVSP_MAIN_GAIA'
506 !
507 
508 ! !==== conversion of thickness to mass ========================!
509  DO icla = 1,nsicla ! convert for each class thickness evolution to mass evolution
510  DO i = 1,npoin
511  evcl_mb%ADR(icla)%P%R(i)
512  & = zfcl_c%ADR(icla)%P%R(i) / mpa2t(icla)
513  mass_sand(icla,1,i) = avail(i,1,icla) / mpa2t(icla)*es(i,1)
514  ENDDO
515  ENDDO
516 
517  ENDIF !VSMTYPE
518 
519  ENDDO ! ENDLOOP ON NOMBLAY
520 !-----------------------------------------------------------------------
521 !
522 100 CONTINUE
523 !
524 !----------------------------------------------------------------------
525 ! TOTAL MASS OF SEDIMENT IN EACH POINT AFTER SLIDING
526 ! JUST FOR CHECK THAT EVOLUTION OF MASS = 0
527 ! EVOLUTIONS FOR EACH CLASS ARE ADDED: TOTAL MASS EVOLUTION
528 ! EVOL_MM%R() CONTAINED -MASS OF SEDIMENT BEFORE SLIDING
529 ! WE ADD MASS OF SEDIMENT AFTER SLIDING TO HAVE EVOLUTION
531  IF(nsand.GE.1) THEN
532  DO ipoin = 1,npoin
533  DO ilayer = 1,nomblay
534  DO isand = 1,nsand
535  evol_mm%R(ipoin) = evol_mm%R(ipoin)
536  & + mass_sand(isand,ilayer,ipoin)
537  ENDDO
538  ENDDO
539  ENDDO
540  ENDIF
541  IF(nmud.GE.1) THEN
542  DO ipoin = 1,npoin
543  DO ilayer = 1,nomblay
544  DO imud = 1,nmud
545  evol_mm%R(ipoin) =evol_mm%R(ipoin)
546  & + mass_mud(imud,ilayer,ipoin)
547  ENDDO
548  ENDDO
549  ENDDO
550  ENDIF
551 !
552  RETURN
553  END
subroutine cvsp_main_gaia(ZFCL_W, ZF, NSICLA, NPOIN)
Definition: cvsp_main_gaia.f:7
subroutine maxslope_gaia
Definition: maxslope_gaia.f:4
type(bief_obj), target zfcl_c
Bed evolution per class (due to bedload)
double precision pi
Pi.
double precision, dimension(:,:,:), allocatable, target mass_mud
Surface mass of mud (kg/m2), for imud,ilayer,ipoin.
type(bief_obj), target evcl_mb
Mass evolution for class (due to bedload)
integer, target nomblay
Number of bed load model layers = NUMSTRAT+1 to take the active layer into account.
double precision, dimension(:), pointer x
2d coordinates of the mesh
type(bief_obj), target zr
Non erodable (rigid) bottom elevation.
type(bief_obj), target volu2d
Integral of bases.
double precision, dimension(:,:), allocatable, target es
Layer thicknesses as double precision.
double precision, dimension(:,:,:), allocatable flux_mass_mud
Flux of mud transfert inbetween layers.
double precision, dimension(:,:,:), allocatable, target mass_sand
Surface mass of sand (kg/m2), for isand,ilayer,ipoin.
integer nsand
Total number of sand.
integer, pointer nelem
Number of elements in the mesh.
type(bief_obj), target tb
Blocks of working arrays.
type(bief_obj), target unsv2d
Inverse of integral of bases.
double precision, target phised
Friction angle of the sediment.
type(bief_obj), target evol_mm
Total mass evolution due to sliding (maxslope) [kg/m2].
integer, target nsicla
Number of sediment classes of bed material (less than NISCLM)
integer vsmtype
For the Continous Vertical Sorting MODEL.
type(bief_obj), pointer t2
integer nmud
Total number of muds.
subroutine cpstvc(X, Y)
Definition: cpstvc.f:7
integer debug
Debugger.
subroutine os(OP, X, Y, Z, C, IOPT, INFINI, ZERO)
Definition: os.f:7
type(bief_obj), target zf
Bottom elevation.
type(bief_mesh), target mesh
Mesh structure.
subroutine parcom(X, ICOM, MESH)
Definition: parcom.f:7
double precision, dimension(:), allocatable, target mpa2t
conversion mass per area to thickness
double precision, dimension(:,:,:), allocatable flux_mass_sand
Flux of sand transfert inbetween layers.
subroutine bed1_update(ZR, ZF, VOLU2D)
Definition: bed1_update.f:7
double precision, dimension(:,:,:), allocatable, target avail
Sediment fraction for each layer, class, point.
integer, pointer npoin
Number of 2d points in the mesh.
type(bief_obj), target v2dpar
Integral of bases in parallel.
Definition: bief.f:3