The TELEMAC-MASCARET system  trunk
lecdon_gaia.f
Go to the documentation of this file.
1 ! **********************
2  SUBROUTINE lecdon_gaia
3 ! **********************
4 !
5  &(motcar,file_desc,path,ncar,code,cas_file,dico_file)
6 !
7 !***********************************************************************
8 ! GAIA
9 !***********************************************************************
10 !
12 !
13 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
21 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
22 !
26 !
28  IMPLICIT NONE
29 !
30 !!-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
31 !
32  INTEGER, INTENT(IN) :: NCAR
33  CHARACTER(LEN=24), INTENT(IN) :: CODE
34  CHARACTER(LEN=PATH_LEN), INTENT(IN) :: PATH
35  CHARACTER(LEN=PATH_LEN), INTENT(INOUT) :: MOTCAR(maxkeyword)
36  CHARACTER(LEN=PATH_LEN), INTENT(INOUT) :: FILE_DESC(4,maxkeyword)
37 ! API
38  CHARACTER(LEN=PATH_LEN), INTENT(IN) :: CAS_FILE
39  CHARACTER(LEN=PATH_LEN), INTENT(IN) :: DICO_FILE
40 !
41 !!-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+-+
42 !
43  INTEGER :: I,K,ERR
44  INTEGER :: MOTINT(maxkeyword)
45  INTEGER :: TROUVE(4,maxkeyword)
46  INTEGER :: ADRESS(4,maxkeyword)
47  INTEGER :: DIMENS(4,maxkeyword)
48  DOUBLE PRECISION :: SUMAVAI
49  DOUBLE PRECISION :: MOTREA(maxkeyword)
50  LOGICAL :: DOC,EFFPEN
51  LOGICAL :: MOTLOG(maxkeyword)
52  CHARACTER(LEN=PATH_LEN) :: NOM_CAS
53  CHARACTER(LEN=PATH_LEN) :: NOM_DIC
54  CHARACTER(LEN=72) :: MOTCLE(4,maxkeyword,2)
55 
56  INTEGER :: ID_DICO, ID_CAS
57 !
58 !-----------------------------------------------------------------------
59 !
60  CHARACTER(LEN=2) CHAR2
61  CHARACTER(LEN=PATH_LEN) TEMPVAR
62 !
63 !-----------------------------------------------------------------------
64 !
65 ! INITIALISES THE VARIABLES FOR DAMOCLES CALL
66 !
67  DO k = 1, maxkeyword
68 ! A FILENAME NOT GIVEN BY DAMOCLES WILL BE RECOGNIZED AS A WHITE SPACE
69 ! (IT MAY BE THAT NOT ALL COMPILERS WILL INITIALISE LIKE THAT)
70  motcar(k)(1:1)=' '
71 !
72  dimens(1,k) = 0
73  dimens(2,k) = 0
74  dimens(3,k) = 0
75  dimens(4,k) = 0
76  ENDDO
77 !
78 ! WRITES OUT INFO
79  doc = .false.
80 !
81 !-----------------------------------------------------------------------
82 ! OPENS DICTIONNARY AND STEERING FILES
83 !-----------------------------------------------------------------------
84 !
85  IF(ncar.GT.0) THEN
86 !
87  nom_dic=path(1:ncar)//'GAIDICO'
88  nom_cas=path(1:ncar)//'GAICAS'
89 !
90  ELSE
91 !
92  nom_dic='GAIDICO'
93  nom_cas='GAICAS'
94 !
95  ENDIF
96  IF((cas_file(1:1).NE.' ').AND.(dico_file(1:1).NE.' ')) THEN
97  WRITE(lu,*) 'FIXED DICO AND STEERING FILE PRESENT'
98  nom_dic=dico_file
99  nom_cas=cas_file
100  WRITE(lu,*) 'NOM_DIC',nom_dic
101  WRITE(lu,*) 'NOM_CAS',nom_cas
102  ENDIF
103 !
104  CALL get_free_id(id_dico)
105  OPEN(id_dico,file=nom_dic,form='FORMATTED',action='READ')
106  CALL get_free_id(id_cas)
107  OPEN(id_cas,file=nom_cas,form='FORMATTED',action='READ')
108 !
109 !-----------------------------------------------------------------------
110 ! CALLS DAMOCLES
111 !-----------------------------------------------------------------------
112 !
113  CALL damocle( adress, dimens ,maxkeyword, doc , lng , lu ,
114  & motint, motrea ,motlog , motcar ,
115  & motcle, trouve ,id_dico, id_cas,.false. ,file_desc)
116 !
117 !-----------------------------------------------------------------------
118 ! CLOSES DICTIONNARY AND STEERING FILES
119 !-----------------------------------------------------------------------
120 !
121  CLOSE(id_dico)
122  CLOSE(id_cas)
123 !
124 ! DECODES 'SUBMIT' CHAINS
125 !
126  CALL read_submit(gai_files,maxlu_gai,file_desc,maxkeyword)
127 !
128 !-----------------------------------------------------------------------
129 !
130 ! RETRIEVES FILES NUMBERS IN GAIA FORTRAN PARAMETERS
131 ! AT THIS LEVEL LOGICAL UNITS ARE EQUAL TO THE FILE NUMBER
132 !
133  DO i=1,maxlu_gai
134  IF(gai_files(i)%TELNAME.EQ.'GAIGEO') THEN
135  gaigeo=i
136  ELSEIF(gai_files(i)%TELNAME.EQ.'GAICLI') THEN
137  gaicli=i
138  ELSEIF(gai_files(i)%TELNAME.EQ.'GAIPRE') THEN
139  gaipre=i
140  ELSEIF(gai_files(i)%TELNAME.EQ.'GAIRES') THEN
141  gaires=i
142  ELSEIF(gai_files(i)%TELNAME.EQ.'GAIREF') THEN
143  gairef=i
144  ELSEIF(gai_files(i)%TELNAME.EQ.'GAICOU') THEN
145  gaicou=i
146  ELSEIF(gai_files(i)%TELNAME.EQ.'GAIFON') THEN
147  gaifon=i
148  ELSEIF(gai_files(i)%TELNAME.EQ.'GAISEC') THEN
149  gaisec=i
150  ELSEIF(gai_files(i)%TELNAME.EQ.'GAISEO') THEN
151  gaiseo=i
152  ELSEIF(gai_files(i)%TELNAME.EQ.'GAILIQ') THEN
153  gailiq=i
154  ELSEIF(gai_files(i)%TELNAME.EQ.'GAIFLX') THEN
155  gaiflx=i
156 ! === NESTOR FILES
157  ELSEIF(gai_files(i)%TELNAME.EQ.'SINACT') THEN
158  sinact=i
159  ELSEIF(gai_files(i)%TELNAME.EQ.'SINPOL') THEN
160  sinpol=i
161  ELSEIF(gai_files(i)%TELNAME.EQ.'SINREF') THEN
162  sinref=i
163  ELSEIF(gai_files(i)%TELNAME.EQ.'SINRST') THEN
164  sinrst=i
165 ! ===
166  ELSEIF(gai_files(i)%TELNAME.EQ.'VSPRES') THEN
167  vspres=i
168 !
169  ENDIF
170  ENDDO
171 !
172 !-----------------------------------------------------------------------
173 !
174 ! ALLOCATING ARRAYS DEPENDING ON MAXFRO
175 !
176 !-----------------------------------------------------------------------
177 !
178 ! MAXIMUM NUMBER OF BOUNDARIES
179  maxfro = motint( adress(1,58) )
180 !
181  ALLOCATE(soldis(maxfro))
182  ALLOCATE(okcgl(maxfro))
183  ALLOCATE(okqgl(maxfro))
184 !
185  DO k=1,maxfro
186  okcgl(k)=.true.
187  okqgl(k)=.true.
188  ENDDO
189 !
190 !-----------------------------------------------------------------------
191 !
192 ! ASSIGNS THE STEERING FILE VALUES TO THE PARAMETER FORTRAN NAME
193 !
194 !-----------------------------------------------------------------------
195 !
196 ! OPTION OF MATRIX ASSEMBLY IS HARD-CODED
197 !
198  optass = 1
199  produc = 1
200 !
201 ! DISCRETISATIONS
202 !
203  ielmt = 11 ! SEDIMENTOLOGICAL VARIABLES
204  ielmh_gai = 11 ! VARIABLES ASSOCIATED WITH WATER DEPTH
205  ielmu_gai = 11 ! VARIABLES ASSOCIATED WITH VELOCITIES
206 !
207 ! FOR NOW PRINTOUTS START AT ZERO
208 !
209  ptinig = 0
210  ptinil = 0
211 !
212 ! NON-EQUILIBIRUM BEDLOAD
213 !
214 ! ICM = MOTINT( ADRESS(1, 1) )
215  icf = motint( adress(1, 2) )
216  bed_model = motint( adress(1, 60) )
217 ! N1 = MOTINT( ADRESS(1, 5) )
218 
219  optban = motint( adress(1, 11) )
220  lvmac = motint( adress(1, 12) )
221  nsous = motint( adress(1, 14) )
222 !
223  mardat(1) = motint( adress(1, 15) )
224  mardat(2) = motint( adress(1, 15) + 1 )
225  mardat(3) = motint( adress(1, 15) + 2 )
226  martim(1) = motint( adress(1, 16) )
227  martim(2) = motint( adress(1, 16) + 1 )
228  martim(3) = motint( adress(1, 16) + 2 )
229 !
230  npriv = motint( adress(1, 23) )
231  nadvar = motint( adress(1,30) )
232 ! NUMBER OF DIRECTIONS FOR DIFFERENTIATED VARIABLES
233  ad_numofdir = motint( adress(1,59) )
234 !
235 ! NCSIZE = MOTINT( ADRESS(1, 24) )
236 ! NUMBER OF PROCESSORS (ALREADY KNOWN FROM THE HYDRODYNAMICS
237 ! MODULE: MUST BE THE SAME, BUT WHEN USING COUPLED MODELS
238 ! IT CAN BE DIFFERENT BY MISTAKE)
239  IF(ncsize.NE.motint(adress(1,24))) THEN
240  WRITE(lu,*) 'DIFFERENT NUMBER OF PARALLEL PROCESSORS:'
241  WRITE(lu,*) 'DECLARED BEFORE (CASE OF COUPLING ?):',ncsize
242  WRITE(lu,*) 'GAIA :',motint(adress(1,24))
243  WRITE(lu,*) 'VALUE ',ncsize,' IS KEPT'
244  ENDIF
245 
246  produc = motint( adress(1, 33) )
247  optass = motint( adress(1, 34) )
248  slopeff = motint( adress(1, 39) )
249  devia = motint( adress(1, 40) )
250  numstrat = motint( adress(1,251) )
251  nsicla = dimens(4,59)
252 !
253 ! INITIALISATION OF ARRAYS WITH DEFAULT VALUES
254 !
255  DO i=1,nsiclm
256  xmvs0(i) = 2650.d0
257  dcla(i) = 0.01d0
258  ac(i) = -9.d0
259  xwc0(i) = -9.d0
260  tocd_mud0(i) = 1000.d0
261  partheniades0(i) = 1.d-3
262  hidi(i) = 1.d0
263  ava0(i) = 0.d0
264  ENDDO
265 ! IF ONLY ONE SEDIMENT AVA0(1) = 1.D0
266  ava0(1) = 1.d0
267 !
268  IF(nsicla.GT.1) THEN
269 ! WHEN THERE IS A SEDIMENT MIXTURE THE
270 ! GRANULOMETRY DISTRIBUTION ALONG THE VERTICAL
271 ! IN THE BED EVOLVES AND THIS IS REPRESENTED
272 ! THROUGH THE USE OF AN ACTIVE LAYER
273 ! NOMBLAY IS THUS THE NUMBER OF PHYSICAL LAYERS
274 ! PLUS ONE, TO ACCOUNT FOR THE ACTIVE LAYER
275  nomblay = numstrat+1
276  ELSE
277  nomblay = numstrat
278  ENDIF
279 !
280  DO i=1,nomblay
281  xkv0(i) = 0.4d0
282  ENDDO
283 !
284 ! FILL DEFAULT LAYER THICKNESSES
285  DO i=1,numstrat
286  sed_thick(i)=100.d0/dble(numstrat)
287  ENDDO
288 !
289 !
290  hidfac = motint( adress(1, 52) )
291  icq = motint( adress(1, 41) )
292 ! CONTROL SECTIONS
293  ncp=dimens(1,42)
294  ALLOCATE(ctrlsc(ncp),stat=err)
295  IF(err.NE.0) THEN
296  WRITE(lu,2039) err
297 2039 FORMAT(1x,'LECDON_GAIA:',/,1x,
298  & 'ERROR DURING ALLOCATION OF CTRLSC: ',/,1x,
299  & 'ERROR CODE: ',1i6)
300  ENDIF
301  IF(ncp.GE.1) THEN
302  DO k=1,ncp
303  ctrlsc(k) = motint( adress(1,42) + k-1 )
304  ENDDO
305  ENDIF
306 ! COORDINATES OF THE ORIGIN
307  i_orig = motint( adress(1,43) )
308  j_orig = motint( adress(1,43)+1 )
309  debug = motint( adress(1,44) )
310 !
311  icr = motint(adress(1,46) )
312 !
313  iks = motint(adress(1,47) )
314 !
315 ! MAXIMUM NUMBER OF ITERATIONS FOR ADVECTION SCHEMES
316 ! TO COMPUTE THE BEDLOAD, FINITE VOLUME SCHEME CENTERED UPWIND, OR
317 ! POSITIVE DEPTH
318  maxadv = motint(adress(1,54))
319 !
320 ! VERTICAL ADVECTION SCHEME
321 !
322  setdep = motint( adress(1,62))
323 !
324  type_houle = motint( adress(1,63) )
325 !
326 ! CVSM VARIABLES
327  vsmtype = motint( adress(1,64))
328  pro_max_max = motint( adress(1,65))
329  cvsmpperiod = motint( adress(1,66))
330  alt_model = motint( adress(1,67))
331 ! As long as there is no coupling period from telemac...
332  percou = 1
333 
334 ! ############### !
335 ! REAL KEYWORDS !
336 ! ############### !
337 !
338  IF(trouve(2,3).GE.1.AND.dimens(2,3).EQ.nsicla) THEN
339  DO k=1,nsicla
340  xmvs0(k) = motrea(adress(2,3)+k-1)
341  ENDDO
342  ELSEIF(trouve(2,3).GE.1.AND.
343  & dimens(2,3).LT.nsicla.AND.dimens(2,3).GT.0) THEN
344 ! READING WHAT HAS BEEN GIVEN
345  DO k=1,dimens(2,3)
346  xmvs0(k) = motrea(adress(2,3)+k-1)
347  ENDDO
348 ! COMPLETING WITH THE LAST GIVEN
349  DO k=dimens(2,3)+1,nsicla
350  xmvs0(k) = motrea(adress(2,3)+dimens(2,3)-1)
351  ENDDO
352  ENDIF
353 !
354  IF(trouve(2,4).GE.1.AND.dimens(2,4).EQ.nsicla) THEN
355  DO k=1,nsicla
356  dcla(k) = motrea(adress(2,4)+k-1)
357  ENDDO
358  ELSEIF(trouve(2,4).GE.1.AND.
359  & dimens(2,4).LT.nsicla.AND.dimens(2,4).GT.0) THEN
360 ! READING WHAT HAS BEEN GIVEN
361  DO k=1,dimens(2,4)
362  dcla(k) = motrea(adress(2,4)+k-1)
363  ENDDO
364 ! COMPLETING WITH THE LAST GIVEN
365  DO k=dimens(2,4)+1,nsicla
366  dcla(k) = motrea(adress(2,4)+dimens(2,4)-1)
367  ENDDO
368  ENDIF
369 !
370 ! SHIELDS NUMBERS
371  IF(trouve(2,6).GE.1.AND.dimens(2,6).EQ.nsicla) THEN
372  DO k=1,nsicla
373  ac(k) = motrea(adress(2,6)+k-1)
374  ENDDO
375  ELSEIF(trouve(2,6).GE.1.AND.
376  & dimens(2,6).LT.nsicla.AND.dimens(2,6).GT.0) THEN
377 ! READING WHAT HAS BEEN GIVEN
378  DO k=1,dimens(2,6)
379  ac(k) = motrea(adress(2,6)+k-1)
380  ENDDO
381 ! COMPLETING WITH THE LAST GIVEN
382  DO k=dimens(2,6)+1,nsicla
383  ac(k) = motrea(adress(2,6)+dimens(2,6)-1)
384  ENDDO
385  ENDIF
386  zero = motrea( adress(2, 9) )
387  vce = motrea( adress(2, 10) )
388  hmin = motrea( adress(2, 11) )
389  beta = motrea( adress(2, 16) )
390 ! SETTLING VELOCITIES (SAME TREATMENT AS SHIELDS NUMBERS)
391 !
392  IF(trouve(2,22).GE.1.AND.dimens(2,22).EQ.nsicla) THEN
393  DO k=1,nsicla
394  xwc0(k) = motrea(adress(2,22)+k-1)
395  ENDDO
396  ELSEIF(trouve(2,22).GE.1.AND.
397  & dimens(2,22).LT.nsicla.AND.dimens(2,22).GT.0) THEN
398 ! READING WHAT HAS BEEN GIVEN
399  DO k=1,dimens(2,22)
400  xwc0(k) = motrea(adress(2,22)+k-1)
401  ENDDO
402 ! COMPLETING WITH THE LAST GIVEN
403  DO k=dimens(2,22)+1,nsicla
404  xwc0(k) = motrea(adress(2,22)+dimens(2,22)-1)
405  ENDDO
406  ENDIF
407 !
408  kspratio = motrea( adress(2, 24) )
409  phised = motrea( adress(2, 25) )
410  beta2 = motrea( adress(2, 26) )
411  bijk = motrea( adress(2, 27) )
412  d90 = motrea( adress(2, 28) )
413 !
414 ! COHESIVE SEDIMENT
415 ! +++++++++++++++++
416 !
417  ncouch_tass = motint( adress(1,45) )
418 !
419 ! DEFAULT VALUES, PREVIOUSLY IN THE DICTIONARY
420 ! EVEN IF THERE ARE NOT 10 LAYERS
421  conc_mud0(1) = 50.d0
422  conc_mud0(2) = 100.d0
423  conc_mud0(3) = 150.d0
424  conc_mud0(4) = 200.d0
425  conc_mud0(5) = 250.d0
426  conc_mud0(6) = 300.d0
427  conc_mud0(7) = 350.d0
428  conc_mud0(8) = 400.d0
429  conc_mud0(9) = 450.d0
430  conc_mud0(10) = 500.d0
431 !
432  IF(dimens(2,32).GT.0) THEN
433  DO k=1,dimens(2,32)
434  conc_mud0(k)=motrea( adress(2,32) + k-1 )
435  ENDDO
436  ENDIF
437 
438 !
439 !
440  IF(dimens(2,34).GT.0) THEN
441  DO k=1,dimens(2,34)
442  toce_mud0(k)=motrea( adress(2,34) + k-1 )
443  ENDDO
444  ENDIF
445 !
446  IF(trouve(2,36).GE.1.AND.dimens(2,36).EQ.nsicla) THEN
447  DO k=1,nsicla
448  tocd_mud0(k) = motrea(adress(2,36)+k-1)
449  ENDDO
450  ELSEIF(trouve(2,36).GE.1.AND.
451  & dimens(2,36).LT.nsicla.AND.dimens(2,36).GT.0) THEN
452 ! READING WHAT HAS BEEN GIVEN
453  DO k=1,dimens(2,36)
454  tocd_mud0(k) = motrea(adress(2,36)+k-1)
455  ENDDO
456 ! COMPLETING WITH THE LAST GIVEN
457  DO k=dimens(2,36)+1,nsicla
458  tocd_mud0(k) = motrea(adress(2,36)+dimens(2,36)-1)
459  ENDDO
460  ENDIF
461 !
462 !
463 ! PARTHENIADES WITH CONVERSION TO M/S
464 !
465 !can be also NOMBLAY
466  IF(trouve(2,37).GE.1.AND.dimens(2,37).EQ.nomblay) THEN
467  DO k=1,nomblay
468  partheniades0(k) = motrea(adress(2,37)+k-1)
469  ENDDO
470  ELSEIF(trouve(2,37).GE.1.AND.
471  & dimens(2,37).LT.nomblay.AND.dimens(2,37).GT.0) THEN
472 ! READING WHAT HAS BEEN GIVEN
473  DO k=1,dimens(2,37)
474  partheniades0(k) = motrea(adress(2,37)+k-1)
475  ENDDO
476 ! COMPLETING WITH THE LAST GIVEN
477  DO k=dimens(2,37)+1,nomblay
478  partheniades0(k) = motrea(adress(2,37)+dimens(2,37)-1)
479  ENDDO
480  ENDIF
481 !
482 ! DEFAULT VALUES, PREVIOUSLY IN THE DICTIONARY
483 ! EVEN IF THERE ARE NOT 10LAYERS
484  trans_mass0(1) = 5.d-5
485  trans_mass0(2) = 4.5d-5
486  trans_mass0(3) = 4.d-5
487  trans_mass0(4) = 3.5d-5
488  trans_mass0(5) = 3.d-5
489  trans_mass0(6) = 2.5d-5
490  trans_mass0(7) = 2.d-5
491  trans_mass0(8) = 1.5d-5
492  trans_mass0(9) = 1.d-5
493  trans_mass0(10) = 0.d0
494 !
495  IF(dimens(2,33).GT.0) THEN
496  DO k=1,dimens(2,33)
497  trans_mass0(k)=motrea( adress(2,33) + k-1 )
498  ENDDO
499  ENDIF
500 !
501 ! BED LAYER THICKNESS
502  IF(trouve(2,42).GE.1.AND.dimens(2,42).EQ.numstrat) THEN
503  DO k=1,numstrat
504  sed_thick(k) = motrea(adress(2,42)+k-1)
505  ENDDO
506  ELSEIF(trouve(2,42).GE.1.AND.
507  & dimens(2,42).LT.numstrat.AND.dimens(2,42).GT.0) THEN
508 ! READING WHAT HAS BEEN GIVEN
509  DO k=1,dimens(2,42)
510  sed_thick(k) = motrea(adress(2,42)+k-1)
511  ENDDO
512 ! COMPLETING WITH THE LAST GIVEN
513  DO k=dimens(2,42)+1,numstrat
514  sed_thick(k) = motrea(adress(2,42)+dimens(2,42)-1)
515  ENDDO
516  WRITE(lu,*)'WARNING:'
517  WRITE(lu,*)'BED LAYER THICKNESS IS NOT SPECIFIED'
518  WRITE(lu,*)'FOR ALL LAYERS'
519  WRITE(lu,*)'LAST VALUE WILL BE KEPT FOR THE'
520  WRITE(lu,*)'LAST LAYERS'
521  ENDIF
522 !
523 ! PRESCRIBED SOLID DISCHARGES
524 !
525  nsoldis=dimens(2,51)
526  IF(nsoldis.GT.0) THEN
527  DO k=1,nsoldis
528  soldis(k)=motrea(adress(2,51)+k-1)
529  ENDDO
530  ENDIF
531 !
532 ! MINIMUM DEPTH FOR BEDLOAD WHEN USING FE
533 !
534  hmin_bedload=motrea(adress(2,52))
535 !
536 ! CENTERING FOR FV
537 !
538  dvf=motrea(adress(2,53))
539 !
540 ! HIDING EXPOSURE MULTI GRAIN MODEL
541 !
542  IF(trouve(2,7).GE.1.AND.dimens(2,7).EQ.nsicla) THEN
543  DO k=1,nsicla
544  hidi(k) = motrea(adress(2,7)+k-1)
545  ENDDO
546  ELSEIF(trouve(2,7).GE.1.AND.
547  & dimens(2,7).LT.nsicla.AND.dimens(2,7).GT.0) THEN
548 ! READING WHAT HAS BEEN GIVEN
549  DO k=1,dimens(2,7)
550  hidi(k) = motrea(adress(2,7)+k-1)
551  ENDDO
552 ! COMPLETING WITH THE LAST GIVEN
553  DO k=dimens(2,7)+1,nsicla
554  hidi(k) = motrea(adress(2,7)+dimens(2,7)-1)
555  ENDDO
556  ENDIF
557 !
558  IF(trouve(2,258).GE.1.AND.dimens(2,258).EQ.nsicla) THEN
559  DO k=1,nsicla
560  ava0(k) = motrea(adress(2,258)+k-1)
561  ENDDO
562  ELSEIF(trouve(2,258).GE.1.AND.
563  & dimens(2,258).LT.nsicla.AND.dimens(2,258).GT.0) THEN
564 ! READING WHAT HAS BEEN GIVEN
565  DO k=1,dimens(2,258)
566  ava0(k) = motrea(adress(2,258)+k-1)
567  ENDDO
568 ! COMPLETING WITH THE LAST GIVEN
569  DO k=dimens(2,258)+1,nsicla
570  ava0(k) = motrea(adress(2,258)+dimens(2,258)-1)
571  ENDDO
572  ENDIF
573 !
574  elay0 = motrea( adress(2,259) )
575 !
576 ! UM: MPM-Factor
577  mpm = motrea( adress(2,260) )
578 ! UM: ALPHA-Factor
579  alpha = motrea( adress(2,261) )
580 ! UM: MOFAC-Factor
581  mofac = motrea( adress(2,262) )
582  mofac_bed = motrea( adress(2,263) )
583 !
584  cgel= motrea( adress(2,29) )
585  cini= motrea( adress(2,35) )
586  turba = motrea(adress(2,40))
587  turbb = motrea(adress(2,41))
588 ! ################## !
589  ! LOGICAL KEYWORDS !
590  ! ################## !
591 ! INDEX 99 IS ALREADY USED FOR KEYWORD 'LIST OF FILES'
592 ! INDEX 54 IS ALREADY USED FOR KEYWORD 'DESCRIPTION OF LIBRARIES'
593 ! INDEX 57 IS ALREADY USED FOR KEYWORD 'DEFAULT EXECUTABLE'
594  ! SPHERICAL EQUATIONS HARD-CODED
595  ! ----------------------------------
596  spheri = .false.
597 !
598  bilma = motlog( adress(3, 1) )
599  IF(bilma) THEN
600  ALLOCATE(sumbedload_b_flux(maxfro))
601  ENDIF
602  bandec = motlog( adress(3, 3) )
603  valid = motlog( adress(3, 4) )
604 ! DTVAR = MOTLOG( ADRESS(3, 5) )
605  susp_sand = motlog( adress(3, 7) )
606  charr = motlog( adress(3, 8) )
607  houle = motlog( adress(3, 10) )
608  const_alayer = motlog( adress(3, 11) )
609  IF(const_alayer.AND.alt_model.NE.0) THEN
610  WRITE(lu,*)'IF CONSTANT ACTIVE LAYER THICKNESS=YES'
611  WRITE(lu,*)'ACTIVE LAYER THICKNESS FORMULA MUST=0'
612  stop
613  ENDIF
614 ! USED TO CHECK GAI_FILES(GAIPRE)%NAME
615  debu = motlog( adress(3, 14) )
616  imp_inflow_c = motlog( adress(3, 15) )
617  seccurrent = motlog( adress(3, 16) )
618  havesecfile = motlog( adress(3, 59) )
619  IF(code(1:9).EQ.'TELEMAC3D') seccurrent = .false.
620  vf = motlog( adress(3, 2) )
621  corr_conv = motlog( adress(3, 18) )
622  slide = motlog( adress(3, 20) )
623  effpen = motlog( adress(3, 22) )
624  IF(.NOT.effpen) THEN
625  slopeff=0
626  devia=0
627  ENDIF
628 !
629 ! COUPLING WITH NESTOR
630  nestor=motlog(adress(3,25))
631 ! V6P1
632  kscalc =motlog(adress(3,26))
633 !
634 ! Settling lag: determines choice between Rouse and Miles concentration profile
635 ! SET_LAG = TRUE : Miles
636 ! = FALSE: Rouse
637 !
638  set_lag = motlog(adress(3,27) )
639 ! Checking the mesh
640  check_mesh = motlog(adress(3,29) )
641 ! NEW IMPLEMENTATION FOR CROSS-SECTION
642  doflux = motlog(adress(3,61) )
643 !##> RK @ BAW
644  IF ( code(1:7) .NE. 'TELEMAC' ) THEN
645 !##> JR @ ADJOINTWARE SYMBOLIC LINEAR SOLVER FOR AD
646 ! SYMBOLIC LINEAR SOLVER FOR AD
647  ad_symblinsolv = motlog( adress(3,30) )
648 ! RESET TANGENTS ON ENTRY IN LINEAR SOLVER CG FOR AD
649  ad_linsolv_resetderiv = motlog( adress(3,31) )
650 ! ALLOW ADDITONAL INTERATIONS IN ITERATIVE LINEAR SOLVERS
651 ! SHOULD NOT BE USED IN PARALLEL MODE
652  ad_linsolv_derivative_convergence = motlog( adress(3,32) )
653 !##< JR @ ADJOINTWARE
654  ENDIF
655 !##< RK @ BAW
656  hinder= motlog( adress(3,9))
657  floc = motlog( adress(3,34))
658 !
659 ! ################################### !
660 ! CHARACTER STRING KEYWORDS !
661 ! ################################### !
662 !
663  titca = motcar( adress(4, 1) )(1:72)
664  sortis = motcar( adress(4, 2) )(1:72)
665  varim = motcar( adress(4, 3) )(1:72)
666  gai_files(gaigeo)%NAME=motcar( adress(4,6) )
667  gai_files(gaicli)%NAME=motcar( adress(4,9) )
668  gai_files(gaipre)%NAME=motcar( adress(4,11) )
669  gai_files(gaires)%NAME=motcar( adress(4,12) )
670  gai_files(gaifon)%NAME=motcar( adress(4,16) )
671  gai_files(gaires)%FMT = motcar( adress(4,31) )(1:8)
672  CALL majus(gai_files(gaires)%FMT)
673 ! RESULT FILE FORMAT FOR PREVIOUS SEDIMENTOLOGICAL
674 ! COMPUTATION...
675  gai_files(gaipre)%FMT = motcar( adress(4,34) )(1:8)
676  CALL majus(gai_files(gaipre)%FMT)
677 ! REFERENCE FILE FORMAT
678  gai_files(gairef)%FMT = motcar( adress(4,33) )(1:8)
679  CALL majus(gai_files(22)%FMT)
680 ! HYDRODYNAMIC FILE FORMAT
681 ! WAVE FILE FORMAT (COUPLING WITH TOMAWAC)
682  gai_files(gaicou)%FMT = motcar( adress(4,35) )(1:8)
683  CALL majus(gai_files(gaicou)%FMT)
684  gai_files(gairef)%NAME=motcar( adress(4,22) )
685 ! NESTOR STEERING FILE
686  gai_files(sinact)%NAME = motcar( adress(4,27) ) ! Nestor
687  gai_files(sinpol)%NAME = motcar( adress(4,40) ) ! Nestor
688  gai_files(sinref)%NAME = motcar( adress(4,41) ) ! Nestor
689  gai_files(sinrst)%NAME = motcar( adress(4,64) ) ! Nestor
690 ! ****** = MOTCAR( ADRESS(4,28) )
691 ! WAVE FILE
692  gai_files(gaicou)%NAME=motcar( adress(4,30) )
693 ! SECTIONS
694  gai_files(gaisec)%NAME=motcar( adress(4,36) )
695  gai_files(gaiseo)%NAME=motcar( adress(4,37) )
696 ! FILE FOR LIQUID BOUNDARIES
697  gai_files(gailiq)%NAME=motcar( adress(4,38) )
698 ! GEOMETRY FILE FORMAT
699  gai_files(gaigeo)%FMT = motcar( adress(4,39) )(1:8)
700  CALL majus(gai_files(gaigeo)%FMT)
701 ! CVSM FILES
702 ! CVSM, But it's not Beautiful
703  tempvar = motcar(adress(4,60) ) !cvsm
705 ! C-VSM RESULTS FILE
706  gai_files(vspres)%NAME=motcar( adress(4,53) ) !cvsm
707  gai_files(vspres)%FMT =motcar( adress(4,55) )(1:8) !cvsm
708  CALL majus(gai_files(vspres)%FMT) !cvsm
709 ! NAMES OF PRIVATE VARIABLES
710  n_names_priv = min(4,dimens(4,42))
711  IF(n_names_priv.GT.0) THEN
712  DO k=1,n_names_priv
713  names_prive(k) = motcar(adress(4,42)+k-1)(1:32)
714  ENDDO
715  ENDIF
716 !
717 ! TYPE OF SEDIMENT AND NSICLA
718  nsand = 0
719  nmud = 0
720  nsusp_tel=0
721  IF(nsicla.GT.0) THEN
722  ALLOCATE(type_sed(nsicla))
723 ! for the moment the following variables are overdimensioned
724  ALLOCATE(num_imud_icla(nsicla))
725  ALLOCATE(num_icla_imud(nsicla))
726  ALLOCATE(num_isand_icla(nsicla))
727  ALLOCATE(num_icla_isand(nsicla))
728  ALLOCATE(num_isusp_icla(nsicla))
729  DO k=1,nsicla
730  num_icla_isand(k)=0
731  num_isand_icla(k)=0
732  num_icla_imud(k)=0
733  num_imud_icla(k)=0
734  num_isusp_icla(k)=0
735  type_sed(k) = motcar(adress(4,59)+k-1)(1:3)
736  IF(type_sed(k).EQ.'CO') THEN
737  sedco(k)=.true.
738  nmud=nmud+1
743 ! Warning: If cohesif sediment forcing suspension
744  susp = .true.
745 !
746  ELSEIF(type_sed(k).EQ.'NCO') THEN
747  sedco(k)=.false.
748  nsand=nsand+1
751  IF(susp_sand.EQV..true.) THEN
754  ENDIF
755  IF(susp_sand.EQV..true.) susp = .true.
756  ELSE
757  WRITE(lu,*)'CHECK TYPE OF SEDIMENT'
758  WRITE(lu,*)'POSSIBLE CHOICES ARE: CO AND NCO'
759  CALL plante(1)
760  stop
761  ENDIF
762  ENDDO
763  ENDIF
764 !
765 ! PRESCRIBED SOLID DISCHARGE DISTRIBUTION
766 !
767  ALLOCATE(ratio_debimp(nsand))
768  nprop = dimens(2,8)
769  IF (nprop.GT.0) THEN
770  DO i = 1, nprop
771  ratio_debimp(i) = motrea(adress(2,8)+i-1)
772  ENDDO
773  ENDIF
774 !
775  IF(nsand.GT.0) THEN
776  IF (.NOT.(charr.OR.susp)) THEN
777  WRITE(lu,*) "YOU HAVE NON COHESIVE SEDIMENT BUT NO "//
778  & "BED LOAD OR SUSPENSION"
779  WRITE(lu,*) "ADD BED LOAD FOR ALL SANDS = YES "
780  WRITE(lu,*) "OR SUSPENSION FOR ALL SANDS = YES "
781  WRITE(lu,*) "IN YOUR STEERING FILE"
782  CALL plante(1)
783  ENDIF
784  ENDIF
785 !
786  IF(nmud.EQ.0) THEN
787  DO i=1,nomblay
788  conc_mud0(i)=0.d0
789  ENDDO
790  ENDIF
791 !
792  IF(nsand.EQ.0.AND.icr.NE.0)THEN
793  WRITE(lu,*)'WITH NO SAND:'
794  WRITE(lu,*)'SKIN FRICTION CORRECTION MUST BE EQUAL TO 0'
795  CALL plante(1)
796  stop
797  ENDIF
798 !
799  IF(numstrat.GE.2.AND.nsand.GT.0) THEN
800  IF(dimens(2,5).NE.nomblay.AND.dimens(2,5).GE.2) THEN ! IF ZERO VALUE of POROSITY : DEFAULT, IF ONE VALUE : FILLS ALL LAYERS
801  WRITE(lu,*)'WARNING:'
802  WRITE(lu,*)'NOT ENOUGH VALUES OF POROSITY GIVEN'
803  WRITE(lu,*)'VALUES WILL BE COMPLETED WITH THE LAST GIVEN'
804  ENDIF
805  ENDIF
806 !
807 ! FILLING THE POROSITY OF STRATUM
808  IF(trouve(2,5).GE.1.AND.dimens(2,5).EQ.nomblay) THEN ! NOMBLAY
809  DO k=1,nomblay
810  xkv0(k) = motrea(adress(2,5)+k-1)
811  ENDDO
812  ELSEIF(trouve(2,5).GE.1.AND.
813  & dimens(2,5).LT.nomblay.AND.dimens(2,5).GT.0) THEN
814 ! READING WHAT HAS BEEN GIVEN
815  DO k=1,dimens(2,5)
816  xkv0(k) = motrea(adress(2,5)+k-1)
817  ENDDO
818 ! COMPLETING WITH THE LAST GIVEN
819  DO k=dimens(2,5)+1,nomblay
820  xkv0(k) = motrea(adress(2,5)+dimens(2,5)-1)
821  ENDDO
822  ENDIF
823 !
824 !
825  n_names_advar = dimens(4,70)
826  nadvar = max( nadvar,n_names_advar ) ! WARNING HERE ?
827  IF(nadvar.GT.0) THEN
828  DO i=1,nadvar
829  WRITE(char2,'(I2)') i
830  names_advar(i) = 'DERIVATIVE '//adjustl(char2)//' '
831  & // '?? '
832  ENDDO
833  ENDIF
834  IF(n_names_advar.GT.0) THEN
835  DO k=1,n_names_advar
836  names_advar(k) = motcar(adress(4,70)+k-1)(1:32)
837  ENDDO
838  ENDIF
839 !
840 ! FLUXLINEFILE
841  gai_files(gaiflx)%NAME=motcar( adress(4,69) )
842 ! C-VSM RESULTS FILE
843  gai_files(vspres)%NAME=motcar( adress(4,53) )
844  gai_files(vspres)%FMT =motcar( adress(4,55) )(1:8)
845  CALL majus(gai_files(vspres)%FMT)
846 !
847  WRITE(lu,102)
848 102 FORMAT(1x,/,19x, '********************************************',/,
849  & 19x, '* LECDON: *',/,
850  & 19x, '* AFTER CALLING DAMOCLES *',/,
851  & 19x, '* CHECKING OF DATA READ *',/,
852  & 19x, '* IN THE STEERING FILE *',/,
853  & 19x, '********************************************',/)
854 !
855 !-----------------------------------------------------------------------
856 !
857 ! LOGICALS FOR OUTPUT VARIABLES
858 !-----------------------------------------------------------------------
859 !
862 !
863  CALL nomvar_gaia
864 !
865  CALL sortie(sortis , mnemo , maxvar , sorleo )
866  CALL sortie(varim , mnemo , maxvar , sorimp )
867 !
868  DO i = 1, 4
869  IF ((npriv.LT.i).AND.(sorleo(i+nvar_priv).OR.
870  & sorimp(i+nvar_priv))) THEN
871  npriv=max(npriv,i)
872  ENDIF
873  ENDDO
874 !
875 !-----------------------------------------------------------------------
876 ! DISABLE OUTPUT OF BEDLOAD SEDIMENT FLOWRATES QS_C, QSCX, QSCY
877 ! IN CASE THERE IS NO BEDLOAD
878 ! THE SEDIMENT FLOWRATES SPECIFIC TO SUSPENSION ARE NOT IN THE
879 ! OUTPUT VARIABLES SO THERE IS NO NEED TO DISABLE THEM
880 !
881  IF(.NOT.charr) THEN
882  sorleo(nvar_qs_c+1)=.false.
883  sorleo(nvar_qsxc+1)=.false.
884  sorleo(nvar_qsyc+1)=.false.
885  sorimp(nvar_qs_c+1)=.false.
886  sorimp(nvar_qsxc+1)=.false.
887  sorimp(nvar_qsyc+1)=.false.
888  ENDIF
889 !
890 !-----------------------------------------------------------------------
891 ! FORCE THE VARIABLES NECESSARY FOR RESTARTING TO BE IN THE OUTPUT
892 ! WITH MUD THE MUD CONCENTRATION, TOCE_MUD, PARTHENIADES AND
893 ! TRANS_MASS (IF BED_MODEL.EQ.2) ARE REQUIRED
894 ! COEFFICIENT IN THE LAYERS
895  IF(nmud.GT.0) THEN
896  DO i=1,nomblay
897  sorleo(nvar_layconc+i)=.true.
898  sorleo(nvar_tocemud+i)=.true.
899  sorleo(nvar_parthe+i)=.true.
900  IF(bed_model.EQ.2) THEN
901  sorleo(nvar_mtrans+i)=.true.
902  ENDIF
903  ENDDO
904  ENDIF
905 !
906 !-----------------------------------------------------------------------
907 !
908 ! INITIALISES MSK (MASKING VARIABLE)
909 ! FOR NOW MASKING IS ONLY DONE FOR ONE 'OPTION FOR THE TREATMENT
910 ! OF TIDAL FLATS'. SHOULD BE OFFERED AS AN OPTION FOR THE USER TO
911 ! CREATE ISLANDS IN THE FUTURE
912  msk = .false.
913  IF (.NOT.bandec) optban=0
914  IF (optban.EQ.2) msk = .true.
915 !
916 !-----------------------------------------------------------------------
917 !
918 ! CHECKS WHETHER THERE IS A REFERENCE FILE
919 !
920  IF (valid.AND.gai_files(gairef)%NAME.EQ.' ') THEN
921  valid=.false.
922  WRITE(lu,71)
923  WRITE(lu,*)
924 71 FORMAT(/,1x,'VALIDATION IS NOT POSSIBLE : ',/
925  & ,1x,'NO REFERENCE FILE ! ')
926  CALL plante(1)
927  stop
928  ENDIF
929 !
930 !MGDL
931 ! CHECKS THE NUMBER OF NSICLA
932 ! IF(NSICLA.GT.10) THEN
933  IF(nsicla.GT.nsiclm) THEN
934  WRITE(lu,81) nsiclm
935  WRITE(lu,*)
936 81 FORMAT(/,1x,'THE MAXIMUM NUMBER OF SEDIMENT CLASSES IS', i2)
937  CALL plante(1)
938  stop
939  ENDIF
940 ! CHECKS THE SUM OF INITIAL AVAI
941  sumavai = 0.d0
942  DO i=1,nsicla
943  sumavai = sumavai + ava0(i)
944  ENDDO
945  IF(abs(sumavai-1).GE.1.d-8) THEN
946  WRITE(lu,91)
947 91 FORMAT(/,1x,'SUM OF SEDIMENT FRACTIONS IS NOT 1 ')
948  WRITE(lu,*)
949  WRITE(lu,*) 'IT IS EQUAL TO', sumavai
950  CALL plante(1)
951  stop
952  ENDIF
953 !
954 ! CHECKS THAT THE BEDLOAD TRANSPORT FORMULATION AND THE HIDING
955 ! FACTOR FORMULATION CAN GO TOGETHER
956 !
957  IF ((hidfac.EQ.3.AND.icf.NE.6).OR.
958  & (hidfac.EQ.1.AND.icf.NE.1).OR.
959  & (hidfac.EQ.2.AND.icf.NE.1)) THEN
960  WRITE(lu,111)
961  WRITE(lu,*)
962 111 FORMAT(/,1x,'CHOICE OF TRANSPORT FORMULA AND HIDING FACTOR FORMU
963  &LATION NOT ALLOWED ')
964  CALL plante(1)
965  stop
966  ENDIF
967 !
968 ! MORPHOLOGICAL FACTOR
969  IF( ((abs(mofac - 1.d0).GT.1e-8) .OR.
970  & (abs(mofac_bed - 1.d0).GT.1e-8))
971  & .AND.bed_model.NE.1)THEN
972  WRITE(lu,*)'MORPHOLOGICAL FACTOR IS NOT RECOMMANDED
973  & WITH BED MODEL AND CONSOLIDATION
974  & ->DISTORSION OF TIME'
975  CALL plante(1)
976  stop
977  ENDIF
978 !
979 ! COMPUTATION CONTINUED
980 !
981  IF(debu) THEN
982  IF(gai_files(gaipre)%NAME(1:1).EQ.' ') THEN
983  WRITE(lu,313)
984 313 FORMAT(/,1x,'COMPUTATION CONTINUED:',/,
985  & 1x,'PREVIOUS SEDIMENTOLOGICAL FILE MISSING')
986  CALL plante(1)
987  stop
988  ENDIF
989  ELSE
990  IF(gai_files(gaipre)%NAME(1:1).NE.' ') THEN
991  gai_files(gaipre)%NAME(1:1)=' '
992  WRITE(lu,213)
993 213 FORMAT(/,1x,'NO COMPUTATION CONTINUED:',/,
994  & 1x,'PREVIOUS SEDIMENTOLOGICAL FILE IGNORED')
995  ENDIF
996  ENDIF
997 !
998 ! METHODS NOT CODED UP FOR SUSPENSION
999 ! -------------------------------------------
1000 !
1001  IF(.NOT.houle) gai_files(gaicou)%NAME(1:1)=' '
1002  IF(houle.AND.charr) THEN
1003  IF(icf.NE.4.AND.icf.NE.5.AND.icf.NE.8.AND.
1004  & icf.NE.9.AND.icf.NE.0) THEN
1005  WRITE(lu,1304) icf
1006 1304 FORMAT(' TRANSPORT FORMULA',1i3,1x,
1007  & 'DOES NOT TAKE WAVES INTO ACCOUNT,',/,1x,
1008  & 'TRY 4, 5, 8 OR 9')
1009  CALL plante(1)
1010  stop
1011  ENDIF
1012  ENDIF
1013 !
1014 ! BEDLOAD AND SUSPENDED TRANSPORT COUPLING
1015 ! ---------------------------------
1016 !
1017  IF((icf==30.OR.icf==3.OR.icf==9).AND.susp.AND.charr) THEN
1018  WRITE(lu,1302) icf
1019  CALL plante(1)
1020  stop
1021  ENDIF
1022 1302 FORMAT('FOR THE FORMULA',1i3,/,1x,
1023  & 'THE SUSPENSION TERM IS CALCULATED TWICE,'
1024  & ,' WITH TOTAL LOAD FORMULA AND SUSPENSION ')
1025 !
1026 !
1027  IF(icq.EQ.2.AND.(.NOT.charr)) THEN
1028  WRITE(lu,1402) icq
1029 1402 FORMAT('FOR THE BIJKER REFERENCE CONCENTRATION',1i3,/,1x,
1030  & 'BEDLOAD MUST BE COMPUTED, CHOOSE:',/,1x,
1031  & 'BEDLOAD = YES')
1032  CALL plante(1)
1033  stop
1034  ENDIF
1035 !
1036 ! CHECKS CONGAITENCY OF BEDLOAD LAWS
1037 !
1038 ! SOULSBY SLOPE EFFECT : REQUIRES A THRESHOLD FORMULA
1039 !
1040  IF(slopeff.EQ.2) THEN
1041 ! CHECK FOR ICF=6
1042 ! IF(ICF.NE.1.AND.ICF.NE.6) THEN
1043  IF(icf.NE.1) THEN
1044  WRITE(lu,1404) icf
1045 1404 FORMAT('BED-LOAD TRANSPORT FORMULA, HERE ICF=',1i3,/,1x,
1046  & 'MUST HAVE A THRESHOLD',/,1x,
1047  & 'IF FORMULA FOR SLOPE EFFECT=2 (SOULSBY)')
1048  ENDIF
1049  ENDIF
1050 !
1051 ! V6P0 : COHERENCE IF CONSOLIDATION MODEL IS USED
1052 ! CSF_VASE STEM FROM THE FIRST LAYER OF THE MULTI-LAYER MODEL
1053 ! +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
1054 ! Mots cles supprimes vitesse critique d'erosion et concentration du lit
1055 ! si premiere couche est vide cela n'est pas correct
1056  IF(nomblay.GT.nlaymax) THEN
1057  WRITE (lu,*) 'NUMBER OF BED LOAD MODEL LAYERS LARGER THAN '
1058  WRITE (lu,*) 'THE MAXIMUM PROGRAMMED VALUE OF ', nlaymax
1059  CALL plante(1)
1060  stop
1061  ENDIF
1062 ! IF(NOMBLAY.LT.2) THEN
1063 ! WRITE (LU,*) 'BEWARE: NUMBER OF BED LOAD MODEL LAYERS'
1064 ! WRITE (LU,*) '======= LOWER THAN THE DEFAULT VALUE OF 2'
1065 ! ENDIF
1066 !
1067 ! ERROR MESSAGE FOR BED MODEL CHOISE
1068 !
1069  IF(bed_model.EQ.3) THEN
1070  WRITE(lu,*)'ONLY BED_MODEL = 1 OR 2 ARE AVAILABLE'
1071  WRITE(lu,*)'BED_MODEL 3: STILL TO DO'
1072  CALL plante(1)
1073  stop
1074  ENDIF
1075 !
1076 !----------------------------------------------------------------
1077 !
1078  RETURN
1079  END
integer iks
Bed roughness predictor option.
integer ncouch_tass
Number of layers for consolidation.
logical bilma
Mass balance.
integer icq
Reference concentration formula.
type(bief_obj), target e
Evolution of the bed mass at each point for each time step [kg/m2].
integer, target gaires
Various files ranks, which are also logical units if no coupling.
double precision hmin
Minimal value of the water height: below this value, the sediment flow rate is set to 0...
logical doflux
Fluxline.
double precision elay0
Wanted active layer thickness; ELAYO is a target value for ELAY, but ELAY may be lower than ELAY0 if ...
character(len=72) titca
Title of the case.
integer produc
Matrix-vector product.
logical seccurrent
Secondary currents.
logical corr_conv
Correction on convection velocity.
logical, target susp
Suspension : yes if there is at least one suspended sediment this is the case if there is mud or if s...
integer, target nomblay
Number of bed load model layers = NUMSTRAT+1 to take the active layer into account.
logical valid
Validation.
subroutine lecdon_split_outputpoints(INT_LIST, POINT_ARRAY, FULLOUTPUT)
integer nsoldis
Number of given solid discharges given by user.
double precision, dimension(:), pointer x
2d coordinates of the mesh
integer, parameter maxlu_gai
Maximum rank of logical units as declared in submit strings in the dictionary.
integer icr
Skin friction correction.
double precision, dimension(nlaymax) sed_thick
Thickness of each bed layer (constant)
type(bief_file), dimension(maxlu_gai), target gai_files
For storing information on files.
integer ptinig
First time from which to write the graphical outputs.
double precision cgel
Weak soil concentration for mud.
logical houle
Include wave effects.
subroutine read_submit(FILES, NFILES, SUBMIT, NMOT)
Definition: read_submit.f:7
integer numstrat
Number of layers of initial stratification.
double precision mofac
Morphological factor on the hydrodynamics: distorts the evolution of the hydrodynamics with respect t...
logical floc
Include floculation effects.
integer bed_model
Bed model (3 choices, cf. dico)
integer, target gaipre
integer, dimension(3) martim
Original hour of time.
double precision, dimension(:), allocatable ratio_debimp
Ratio of sand in the prescribed solid discharge.
double precision cini
Threshold concentration for hindered settling.
integer i_orig
Coordinates of the origin.
logical, dimension(maxvar) sorimp
Listing output.
integer optban
Option for the treatment of tidal flats.
character(len=32), dimension(maxvar) names_advar
Names of differenting arrays (given by user)
logical, dimension(maxvar) sorleo
Graphical output.
integer, parameter maxkeyword
double precision vce
Water viscosity: it is defined here because the viscosity set in TELEMAC2D or TELEMAC3D may not b the...
double precision turbb
Coefficient for floc destruction.
integer, dimension(:), allocatable num_isusp_icla
Tables to switch from suspended sediment number to class number.
logical cvsm_out_full
C-VSM_FULL WRITES OUT (OR NOT) EVER.
integer, parameter maxvar
Maximum number of output variables.
double precision, dimension(nlaymax) toce_mud0
Critical erosion shear stress of the mud per layer.
integer optass
Matrix storage.
integer, dimension(100) cvsmoutput
CHOOSE POINTS or FULL MODEL AS PRINTOUT.
double precision, dimension(:), allocatable soldis
Prescribed solid discharges.
integer slopeff
Formula for slope effect.
double precision, target alpha
Secondary Current Alpha Coefficient.
double precision, dimension(nsiclm) xmvs0
Sand density.
double precision, dimension(nsiclm) tocd_mud0
Critical shear stress for mud deposition.
integer nprop
Number of class proportion for imposed discharge given by user.
integer setdep
Options for the vertical advection-diffusion scheme with settling velocity (only relevant in 3D)...
logical imp_inflow_c
Imposed concentration in inflow.
integer nsand
Total number of sand.
logical slide
Sediment slide.
double precision zero
Parameter used for clipping variables or testing values against zero.
integer type_houle
Type of waves (regular or irregular)
logical havesecfile
Secondary currents radii file.
subroutine lecdon_gaia(MOTCAR, FILE_DESC, PATH, NCAR, CODE, CAS_FILE, DICO_FILE)
Definition: lecdon_gaia.f:7
integer, target gaicou
integer, dimension(:), allocatable num_isand_icla
integer icf
Bed-load transport formula.
integer ptinil
First time from which to write the listing outputs.
double precision, dimension(nlaymax) partheniades0
Partheniades erosion coefficient: depends on the type of erosion so it actually varies on the sedimen...
double precision dvf
Upwinding for Exner FV.
integer alt_model
CHOOSE A MODEL FOR ESTIMATION OF A DYNAMIC ACTIVE LAYER THICKNESS.
double precision turba
Flocculation coefficient.
integer, dimension(:), allocatable num_imud_icla
Tables to switch from mud number to class number and from sand number to class number.
character(len=3), dimension(:), allocatable type_sed
Type of sediment (co or nco)
integer, parameter nsiclm
Maximum number of sediment classes.
double precision, dimension(nsiclm), target dcla
Sediment diameter for each class It is only relevant for non-cohesive sediments. For the bedload...
double precision mofac_bed
Morphological factor on the bed: distorts the evolution of the morphodynamics with respect to the hyd...
double precision, dimension(nsiclm) ava0
Initial fraction of each sediment class, the sum of AVA0 over all classes has to be equal to 1...
logical, dimension(:), allocatable okqgl
Used in function qgl_gaia!
logical, dimension(nsiclm) sedco
Cohesive sediments (for each class)
double precision, target phised
Friction angle of the sediment.
logical vf
Use a Finite Volumes formulation.
integer, target nsicla
Number of sediment classes of bed material (less than NISCLM)
logical spheri
Work in spherical coordinates (hard-coded)
integer vsmtype
For the Continous Vertical Sorting MODEL.
logical nestor
Coupling with NESTOR.
subroutine sortie(CHAINE, MNEMO, NBRE, SORLEO)
Definition: sortie.f:7
integer cvsmpperiod
Printout Period for Full Vertical Sorting Model: PRO_D & PRO_F.
double precision, dimension(nlaymax), target xkv0
Initial porosity by layers.
integer lvmac
Vector length (for vectorisation) from steering file.
logical msk
Include masking.
integer nsous
Number of sub-iterations.
logical hinder
Hindered settling switch.
logical bandec
Include tidal flats in the simulation.
logical const_alayer
Constant active layer thickness.
integer nmud
Total number of muds.
integer, parameter nlaymax
Maximum number of layers on the mesh.
subroutine damocle(ADRESS, DIMENS, NMAX, DOC, LLNG, LLU, MOTINT, MOTREA, MOTLOG, MOTCAR, MOTCLE, TROUVE, NFICMO, NFICDA, GESTD, MOTATT)
Definition: damocle.f:9
logical, target susp_sand
Suspension for all sands (mud is assumed to be suspended)
logical debu
Computation continued.
integer, target gaigeo
integer, dimension(:), allocatable num_icla_imud
Tables to switch from class number to mud number and from class number to sand number.
integer maxfro
Maximum number of (liquid boundaries, solid boundaries)
double precision, dimension(nlaymax) conc_mud0
Mud concentration for each bed layer (constant)
character(len=72) sortis
List of the variable to ouput in the result file.
double precision, dimension(nsiclm), target ac
Critical shields parameter.
double precision, target d90
Sediment diameter D90, for sand when only.
integer devia
Formula for deviation.
integer ad_numofdir
Number of directions for differentiating in vector modes.
integer maxadv
Maximum number of iterations used for ensuring positive layer thickness. Beware that with larger time...
double precision bijk
B value for the Bijker formula.
integer percou
COUPLING PERIOD USED IN CVSM TO CALCULATE THE TIME, SHOULD COME FROM TELEMAC ONE DAY.
integer debug
Debugger.
integer pro_max_max
Maximum Number of Profile SECTIONS.
double precision, target kspratio
Ratio between skin friction and mean diameter.
character(len=32), dimension(4) names_prive
Names of private arrays (given by user)
double precision, dimension(nlaymax) trans_mass0
Mass transfer for consolidation between layers.
logical, dimension(:), allocatable okcgl
Used in function cgl_gaia.
integer, dimension(:), allocatable num_icla_isand
double precision beta
Beta coefficient for Koch and Flokstra slope effect formulation.
logical set_lag
Settling lag: determines choice between Rouse and Miles concentration profile SET_LAG = TRUE : Miles ...
integer ielmt
Missing comment.
logical, target charr
Include bedload in the simulation.
logical kscalc
Bed roughness prediction.
character(len=8), dimension(maxvar) mnemo
Mnemo of variables for graphic printouts (b for bottom, etc.)
integer ncp
Number of control sections points.
integer, target gairef
subroutine nomvar_gaia
Definition: nomvar_gaia.f:4
double precision, dimension(nsiclm) hidi
Hiding factor for each sediment class Used only if HIDFAC is set to 0. By default it is set to 1...
subroutine majus(CHAINE)
Definition: majus.f:7
integer, target gaicli
double precision, dimension(:), allocatable sumbedload_b_flux
Sum over classes of bedload boundary flux or cumulated bedload:
character(len=72) varim
List of the variable to print to the listing.
integer npriv
Number of private arrays, number of private arrays with given name.
double precision, target beta2
Parameter for deviation.
double precision hmin_bedload
Minimum depth for bedload.
double precision, target mpm
Meyer Peter Mueller-Coefficient.
integer nadvar
Number of differentiating arrays, and those with a given name.
integer nsusp_tel
Number of suspension sediment classes for TELEMAC3D or TELEMAC2D (less than NISCLM) ...
double precision, dimension(nsiclm) xwc0
Initial settling velocities.
integer, dimension(:), allocatable ctrlsc
Array containing the global number of the points in the control sections.
integer hidfac
Hiding factor formulas.
integer, dimension(3) mardat
Orginal date of time.