The TELEMAC-MASCARET system  trunk
lecdon_tomawac.f
Go to the documentation of this file.
1 ! *************************
2  SUBROUTINE lecdon_tomawac
3 ! *************************
4 !
5  &(file_desc,path,ncar,cas_file,dico_file)
6 !
7 !***********************************************************************
8 ! TOMAWAC V7P1
9 !***********************************************************************
10 !
11 !brief READS THE STEERING FILE THROUGH A DAMOCLES CALL.
12 !
13 !history G.MATTAROLO (EDF)
14 !+ 16/05/2011
15 !+ V6P1
16 !+ Declaration of new keywords defined by
17 !+ E. GAGNAIRE-RENOU for solving new source terms models.
18 !
19 !history G.MATTAROLO (EDF - LNHE)
20 !+ 20/06/2011
21 !+ V6P1
22 !+ Translation of French names of the variables in argument
23 !
24 !history G.MATTAROLO (EDF - LNHE)
25 !+ 25/06/2012
26 !+ V6P2
27 !+ Declaration of new keywords for representing diffraction
28 !
29 !history J-M HERVOUET (EDF R&D LNHE)
30 !+ 01/02/2013
31 !+ V6P3
32 !+ New keywords added. Call to tomawac_constants added.
33 !
34 !history J-M HERVOUET (EDF R&D LNHE)
35 !+ 09/05/2014
36 !+ V7P0
37 !+ Retrieving MODASS for new parallel assembly with integers.
38 !
39 !history J-M HERVOUET (EDF R&D LNHE)
40 !+ 12/09/2014
41 !+ V7P0
42 !+ Retrieving VEGETATION for VEGETATION TAKEN INTO ACCOUNT
43 !
44 !history J-M HERVOUET (EDF LAB, LNHE)
45 !+ 18/05/2015
46 !+ V7P1
47 !+ Adding CHECK_MESH for the keyword 'CHECKING THE MESH'
48 !
49 !history Y AUDOUIN (LNHE)
50 !+ 25/05/2015
51 !+ V7P0
52 !+ Modification to comply with the hermes module
53 !
54 !history J-M HERVOUET (EDF LAB, LNHE)
55 !+ 08/03/2016
56 !+ V7P2
57 !+ Retrieving NPLEO in a different way.
58 !
59 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
60 !| FILE_DESC |-->| STORES THE FILES 'SUBMIT' ATTRIBUTES
61 !| | | IN DICTIONARIES. IT IS FILLED BY DAMOCLES.
62 !| NCAR |-->| LENGTH OF PATH
63 !| PATH |-->| NAME OF CURRENT DIRECTORY
64 !~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
65 !
66  USE bief
70 !
71  USE interface_tomawac, ex_lecdon_tomawac => lecdon_tomawac
72  IMPLICIT NONE
73 !
74 !
75 !-----------------------------------------------------------------------
76 !
77  CHARACTER(LEN=PATH_LEN), INTENT(INOUT) :: FILE_DESC(4,maxkeyword)
78  INTEGER, INTENT(IN) :: NCAR
79  CHARACTER(LEN=PATH_LEN), INTENT(IN) :: PATH
80  CHARACTER(LEN=PATH_LEN), INTENT(IN) :: DICO_FILE
81  CHARACTER(LEN=PATH_LEN), INTENT(IN) :: CAS_FILE
82 !
83 !-----------------------------------------------------------------------
84 !
85 ! ARRAYS USED IN THE DAMOCLES CALL
86 !
87  CHARACTER(LEN=8) MNEMO(maxvar)
88  INTEGER K
89  INTEGER ADRESS(4,maxkeyword),DIMEN(4,maxkeyword)
90  DOUBLE PRECISION MOTREA(maxkeyword)
91  INTEGER MOTINT(maxkeyword)
92  LOGICAL MOTLOG(maxkeyword)
93  CHARACTER(LEN=PATH_LEN) MOTCAR(maxkeyword)
94  CHARACTER(LEN=72) MOTCLE(4,maxkeyword,2)
95  INTEGER TROUVE(4,maxkeyword)
96  LOGICAL DOC
97  CHARACTER(LEN=PATH_LEN) :: NOM_CAS
98  CHARACTER(LEN=PATH_LEN) :: NOM_DIC
99 ! ARGUMENTS
100  INTEGER :: I
101  INTEGER :: ID_DICO,ID_CAS
102 !
103 ! END OF DECLARATIONS FOR DAMOCLES CALL
104 !
105 !
106 !***********************************************************************
107 !
108  WRITE(lu,2)
109 2 FORMAT(1x,/,19x, '********************************************',/,
110  & 19x, '* SUBROUTINE LECDON_TOMAWAC *',/,
111  & 19x, '* CALL OF DAMOCLES *',/,
112  & 19x, '* VERIFICATION OF READ DATA *',/,
113  & 19x, '* ON STEERING FILE *',/,
114  & 19x, '********************************************',/)
115 !
116 !-----------------------------------------------------------------------
117 !
118 ! INITIALISES THE VARIABLES FOR DAMOCLES CALL :
119 !
120  DO k=1,maxkeyword
121 ! A FILENAME NOT GIVEN BY DAMOCLES WILL BE RECOGNIZED AS A WHITE SPACE
122 ! (IT MAY BE THAT NOT ALL COMPILERS WILL INITIALISE LIKE THAT)
123  motcar(k)(1:1)=' '
124 !
125  dimen(1,k) = 0
126  dimen(2,k) = 0
127  dimen(3,k) = 0
128  dimen(4,k) = 0
129  ENDDO
130 !
131 ! WRITES OUT INFO
132  doc = .false.
133 !
134 !-----------------------------------------------------------------------
135 ! OPENS DICTIONNARY AND STEERING FILES
136 !-----------------------------------------------------------------------
137 !
138  IF(ncar.GT.0) THEN
139 !
140  nom_dic=path(1:ncar)//'WACDICO'
141  nom_cas=path(1:ncar)//'WACCAS'
142 !
143  ELSE
144 !
145  nom_dic='WACDICO'
146  nom_cas='WACCAS'
147 !
148  ENDIF
149  IF((cas_file(1:1).NE.' ').AND.(dico_file(1:1).NE.' ')) THEN
150  nom_dic=dico_file
151  nom_cas=cas_file
152  ENDIF
153 !
154  CALL get_free_id(id_dico)
155  OPEN(id_dico,file=nom_dic,form='FORMATTED',action='READ')
156  CALL get_free_id(id_cas)
157  OPEN(id_cas,file=nom_cas,form='FORMATTED',action='READ')
158 !
159  CALL damocle
160  &( adress, dimen , maxkeyword , doc , lng , lu , motint,
161  & motrea, motlog, motcar, motcle , trouve, id_dico, id_cas,
162  & .false.,file_desc)
163 
164  CLOSE(id_dico)
165  CLOSE(id_cas)
166 !
167 ! DECODES 'SUBMIT' CHAINS
168 !
169  CALL read_submit(wac_files,maxlu_wac,file_desc,maxkeyword)
170 !
171 !-----------------------------------------------------------------------
172 !
173 ! RETRIEVES FILE NUMBERS FROM TOMAWAC FORTRAN PARAMETERS
174 !
175  DO i=1,maxlu_wac
176  IF(wac_files(i)%TELNAME.EQ.'WACGEO') THEN
177  wacgeo=i
178  lugeo=>wac_files(i)%LU
179  fmtgeo=> wac_files(i)%FMT
180  ELSEIF(wac_files(i)%TELNAME.EQ.'WACCAS') THEN
181  waccas=i
182  ELSEIF(wac_files(i)%TELNAME.EQ.'WACCLI') THEN
183  waccli=i
184  ELSEIF(wac_files(i)%TELNAME.EQ.'WACFON') THEN
185  wacfon=i
186  lufon=>wac_files(i)%LU
187  namfon=> wac_files(i)%NAME
188  ELSEIF(wac_files(i)%TELNAME.EQ.'WACRES') THEN
189  wacres=i
190  lures=>wac_files(i)%LU
191  namres=> wac_files(i)%NAME
192  fmtres=> wac_files(i)%FMT
193  ELSEIF(wac_files(i)%TELNAME.EQ.'WACREF') THEN
194  wacref=i
195  luref=>wac_files(i)%LU
196  fmtref=> wac_files(i)%FMT
197  ELSEIF(wac_files(i)%TELNAME.EQ.'WACSPE') THEN
198  wacspe=i
199  luspe=>wac_files(i)%LU
200  namspe=> wac_files(i)%NAME
201  ELSEIF(wac_files(i)%TELNAME.EQ.'WACLEO') THEN
202  wacleo=i
203  luleo=>wac_files(i)%LU
204  namleo=> wac_files(i)%NAME
205  fmtleo=> wac_files(i)%FMT
206  ELSEIF(wac_files(i)%TELNAME.EQ.'WACPRE') THEN
207  wacpre=i
208  lupre=>wac_files(i)%LU
209  fmtpre=> wac_files(i)%FMT
210  ELSEIF(wac_files(i)%TELNAME.EQ.'WACRBI') THEN
211  wacrbi=i
212  lurbi=>wac_files(i)%LU
213  fmtrbi=> wac_files(i)%FMT
214  ELSEIF(wac_files(i)%TELNAME.EQ.'WACCOB') THEN
215  waccob=i
216  lucob=>wac_files(i)%LU
217  namcob=> wac_files(i)%NAME
218  fmtcob=> wac_files(i)%FMT
219  ELSEIF(wac_files(i)%TELNAME.EQ.'WACCOF') THEN
220  waccof=i
221  lucof=>wac_files(i)%LU
222  namcof=> wac_files(i)%NAME
223  fmtcof=> wac_files(i)%FMT
224  ELSEIF(wac_files(i)%TELNAME.EQ.'WACBI1') THEN
225  wacbi1=i
226  lubi1=>wac_files(i)%LU
227  fmtbi1=> wac_files(i)%FMT
228  ELSEIF(wac_files(i)%TELNAME.EQ.'WACFO1') THEN
229  wacfo1=i
230  lufo1=>wac_files(i)%LU
231  ELSEIF(wac_files(i)%TELNAME.EQ.'WACVEB') THEN
232  wacveb=i
233  luveb=>wac_files(i)%LU
234  namveb=> wac_files(i)%NAME
235  fmtveb=> wac_files(i)%FMT
236  ELSEIF(wac_files(i)%TELNAME.EQ.'WACVEF') THEN
237  wacvef=i
238  luvef=>wac_files(i)%LU
239  namvef=> wac_files(i)%NAME
240  fmtvef=> wac_files(i)%FMT
241  ELSEIF(wac_files(i)%TELNAME.EQ.'WACMAB') THEN
242  wacmab=i
243  lumab=>wac_files(i)%LU
244  nammab=> wac_files(i)%NAME
245  fmtmab=> wac_files(i)%FMT
246  ELSEIF(wac_files(i)%TELNAME.EQ.'WACMAF') THEN
247  wacmaf=i
248  lumaf=>wac_files(i)%LU
249  nammaf=> wac_files(i)%NAME
250  fmtmaf=> wac_files(i)%FMT
251  ELSEIF(wac_files(i)%TELNAME.EQ.'LEOWXY') THEN
252  leowxy=i
253  namwxy=> wac_files(i)%NAME
254  luwxy=>wac_files(i)%LU
255  ELSEIF(wac_files(i)%TELNAME.EQ.'LEOIXY') THEN
256  leoixy=i
257  namixy=> wac_files(i)%NAME
258  luixy=>wac_files(i)%LU
259  ELSEIF(wac_files(i)%TELNAME.EQ.'IMPSPE') THEN
260  impspe=i
261  ENDIF
262  ENDDO
263 !
264 !-----------------------------------------------------------------------
265 !
266 ! SETTING CONSTANTS (PI, GRAVITY, ETC.)
267 !
268  CALL tomawac_constants
269 !
270 !-----------------------------------------------------------------------
271 !
272 ! ASSIGNS THE STEERING FILE VALUES TO THE PARAMETER FORTRAN NAME
273 !
274 !-----------------------------------------------------------------------
275 !
276 ! INTEGER KEYWORDS
277 !
278  graprd = motint( adress(1, 1) )
279  lisprd = motint( adress(1, 2) )
280  nit = motint( adress(1, 3) )
281  ndire = motint( adress(1, 4) )
282  nf = motint( adress(1, 5) )
283  gradeb = motint( adress(1, 6) )
284  lisfon = motint( adress(1, 7) )
285  svent = motint( adress(1, 8) )
286  smout = motint( adress(1, 9) )
287  sfrot = motint( adress(1, 10) )
288  strif = motint( adress(1, 11) )
289  indic = motint( adress(1, 12) )
290  indiv = motint( adress(1, 13) )
291  nsits = motint( adress(1, 14) )
292  inispe = motint( adress(1, 15) )
293 ! DISSIPATION BY STRONG CURRENT
294  sdscu = motint( adress(1, 16) )
295  nptt = motint( adress(1, 17) )
296  lvmac = motint( adress(1, 18) )
297  sbrek = motint( adress(1, 19) )
298  iqbbj = motint( adress(1, 20) )
299  ihmbj = motint( adress(1, 21) )
300  ifrbj = motint( adress(1, 22) )
301  iwhtg = motint( adress(1, 23) )
302  ifrtg = motint( adress(1, 24) )
303  idisro = motint( adress(1, 25) )
304  iexpro = motint( adress(1, 26) )
305  ifrro = motint( adress(1, 27) )
306  ifrih = motint( adress(1, 28) )
307  ndtbrk = motint( adress(1, 29) )
308  limit = motint( adress(1, 30) )
309 !GM V6P1 - NEW SOURCE TERMS
310  lvent = motint( adress(1, 31) )
311 !GM Fin
312  stria = motint( adress(1, 32) )
313  limspe = motint( adress(1, 33) )
314  lam = motint( adress(1, 34) )
315  indim = motint( adress(1, 35) )
316  idhma = motint( adress(1, 36) )
317  frabi = motint( adress(1, 37) )
318  npriv = motint( adress(1, 38) )
319  frabl = motint( adress(1, 39) )
320 ! COORDINATES OF THE ORIGIN IN (X, Y)
321  i_orig = motint( adress(1, 40) )
322  j_orig = motint( adress(1, 40)+1 )
323 ! DEBUG KEYWORD
324  debug = motint( adress(1, 41) )
325 ! GQM PARAMETERS
326  iq_om1 = motint( adress(1, 42) )
327  nq_te1 = motint( adress(1, 43) )
328  nq_om2 = motint( adress(1, 44) )
329 !
330 ! Diffraction
331  diffra = motint( adress(1, 45) )
332  nptdif = motint( adress(1, 46) )
333 ! END diffraction
334 !
335  diaghf = motint( adress(1, 47) )
336 !
337 ! OPTION FOR SECOND DERIVATIVES
338 !
339  optder = motint( adress(1, 48) )
340 !
341 ! 49 IS PARALLEL PROCESSORS
342 !
343 ! PARALLEL ASSEMBLY MODE
344 !
345  modass = motint( adress(1, 50) )
346 !
347 ! FOR BAJ MODELISATION
348 ! SEE NEXT BAJ FOR CONSEQUENCES ON PARAMETERS
349  cbaj = motint( adress(1, 51) )
350 !
351 ! REAL KEYWORDS
352 !
353  dt = motrea( adress(2, 1) )
354  f1 = motrea( adress(2, 2) )
355  raisf = motrea( adress(2, 3) )
356  IF(dimen(2,4).NE.dimen(2,5)) THEN
357  WRITE(lu,*) 'ABSCISSAE AND ORDINATES OF SPECTRUM PRINTOUT'
358  WRITE(lu,*) 'POINTS MUST BE GIVEN IN EQUAL NUMBERS'
359  WRITE(lu,*) 'THERE ARE HERE',dimen(2,4),' ABCISSAE AND '
360  WRITE(lu,*) dimen(2,5),' ORDINATES'
361  CALL plante(1)
362  stop
363  ENDIF
364  IF(trouve(2,4).EQ.2) THEN
365  npleo = dimen(2,4)
366  ELSE
367 ! IF KEYWORD NOT FOUND
368  npleo = 0
369  ENDIF
370  ALLOCATE(noleo(npleo))
371  ALLOCATE(xleo(npleo))
372  ALLOCATE(yleo(npleo))
373  DO k=1,dimen(2,4)
374  xleo(k)= motrea( adress(2, 4) + k-1)
375  ENDDO
376  DO k=1,dimen(2,5)
377  yleo(k)= motrea( adress(2, 5) + k-1)
378  ENDDO
379  ddc = motrea( adress(2, 6) )
380  cfrot1 = motrea( adress(2, 7) )
381  cmout1 = motrea( adress(2, 8) )
382  cmout2 = motrea( adress(2, 9) )
383  roair = motrea( adress(2, 10) )
384  roeau = motrea( adress(2, 11) )
385  betam = motrea( adress(2, 12) )
386  xkappa = motrea( adress(2, 13) )
387  alpha = motrea( adress(2, 14) )
388  decal = motrea( adress(2, 15) )
389  zvent = motrea( adress(2, 16) )
390  cdrag = motrea( adress(2, 17) )
391  hm0 = motrea( adress(2, 18) )
392  fpic = motrea( adress(2, 19) )
393  gamma = motrea( adress(2, 20) )
394  sigmaa = motrea( adress(2, 21) )
395  sigmab = motrea( adress(2, 22) )
396  alphil = motrea( adress(2, 23) )
397  fetch = motrea( adress(2, 24) )
398  fremax = motrea( adress(2, 25) )
399  teta1 = motrea( adress(2, 26) )*degrad
400  spred1 = motrea( adress(2, 27) )
401  teta2 = motrea( adress(2, 28) )*degrad
402  spred2 = motrea( adress(2, 29) )
403  xlamda = motrea( adress(2, 30) )
404  tailf = motrea( adress(2, 31) )
405  e2fmin = motrea( adress(2, 32) )
406  alfabj = motrea( adress(2, 33) )
407  gambj1 = motrea( adress(2, 34) )
408  gambj2 = motrea( adress(2, 35) )
409  boretg = motrea( adress(2, 36) )
410  gamatg = motrea( adress(2, 37) )
411  alfaro = motrea( adress(2, 38) )
412  gamaro = motrea( adress(2, 39) )
413  gam2ro = motrea( adress(2, 40) )
414  betaih = motrea( adress(2, 41) )
415  em2sih = motrea( adress(2, 42) )
416  coefhs = motrea( adress(2, 43) )
417  xdtbrk = motrea( adress(2, 44) )
418  xlamd = motrea( adress(2, 45) )
419  zrepos = motrea( adress(2, 46) )
420  alflta = motrea( adress(2, 47) )
421  rfmlta = motrea( adress(2, 48) )
422  kspb = motrea( adress(2, 49) )
423  bdispb = motrea( adress(2, 50) )*degrad
424  bdsspb = motrea( adress(2, 51) )*degrad
425  hm0l = motrea( adress(2, 52) )
426  fpicl = motrea( adress(2, 53) )
427  sigmal = motrea( adress(2, 54) )
428  sigmbl = motrea( adress(2, 55) )
429  aphill = motrea( adress(2, 56) )
430  fetchl = motrea( adress(2, 57) )
431  fpmaxl = motrea( adress(2, 58) )
432  teta1l = motrea( adress(2, 59) )*degrad
433  spre1l = motrea( adress(2, 60) )
434  teta2l = motrea( adress(2, 61) )*degrad
435  spre2l = motrea( adress(2, 62) )
436  xlamdl = motrea( adress(2, 63) )
437  gammal = motrea( adress(2, 64) )
438  promin = motrea( adress(2, 65) )
439  vx_cte = motrea( adress(2, 66) )
440  vy_cte = motrea( adress(2, 67) )
441  cimpli = motrea( adress(2, 68) )
442  coefwd = motrea( adress(2, 69) )
443  coefwe = motrea( adress(2, 70) )
444  coefwf = motrea( adress(2, 71) )
445  coefwh = motrea( adress(2, 72) )
446  cmout3 = motrea( adress(2, 73) )
447  cmout4 = motrea( adress(2, 74) )
448  cmout5 = motrea( adress(2, 75) )
449  cmout6 = motrea( adress(2, 76) )
450  seuil = motrea( adress(2, 77) )
451  seuil1 = motrea( adress(2, 78) )
452  seuil2 = motrea( adress(2, 79) )
453  f2difm = motrea( adress(2, 80) )
454 ! TIME UNITS IN FILES
455  unitcob= motrea( adress(2, 81) )
456  unitmab= motrea( adress(2, 82) )
457  unitveb= motrea( adress(2, 83) )
458  unitspe= motrea( adress(2, 88) )
459 ! TIME SHIFTS IN FILES
460  phascob= motrea( adress(2, 84) )
461  phasmab= motrea( adress(2, 85) )
462  phasveb= motrea( adress(2, 86) )
463  phasspe= motrea( adress(2, 89) )
464 ! DISSIPATION COEFFICIENT FOR STRONG CURRENT
465  cdscur = motrea( adress(2, 87) )
466  smax = motrea( adress(2, 90) )
467 !
468 ! LOGICAL KEYWORDS
469 !
470  tsou = motlog( adress(3, 1) )
471  sphe = motlog( adress(3, 2) )
472  proinf = motlog( adress(3, 5) )
473  cousta = motlog( adress(3, 6) )
474  vent = motlog( adress(3, 7) )
475  dontel = motlog( adress(3, 8) )
476  prop = motlog( adress(3, 9) )
477  vensta = motlog( adress(3, 10) )
478  valid = motlog( adress(3, 11) )
479  maree = motlog( adress(3, 12) )
480  trigo = motlog( adress(3, 13) )
481  speuli = motlog( adress(3, 14) )
482  fltdif = motlog( adress(3, 15) )
483  raztim = motlog( adress(3, 16) )
484  vegetation = motlog( adress(3, 17) )
485  check_mesh = motlog( adress(3, 18) )
486  source_on_bnd = motlog( adress(3, 19) )
487  avant = motlog( adress(3, 20) )
488  ecret = motlog( adress(3, 21) )
489  porous = motlog( adress(3, 22) )
490  partel_concat = motlog( adress(3,23) )
491  IF(ncsize.LE.1) partel_concat=.false.
492 !
493 ! STRING KEYWORDS
494 !
495  titcas = motcar( adress(4, 1) ) (1:72)
496  sort2d = motcar( adress(4, 2) ) (1:72)
497 !
498 ! FILES IN THE STEERING FILE
499 !
500  wac_files(wacgeo)%NAME=motcar( adress(4,3) )
501  IF(wac_files(wacgeo)%NAME(1:1).EQ.' ') THEN
502  WRITE(lu,*) 'THE FOLLOWING KEYWORD IS MANDATORY:'
503  WRITE(lu,*) 'GEOMETRY FILE (FICHIER DE GEOMETRIE)'
504  CALL plante(1)
505  stop
506  ENDIF
507 ! NOMFOR = MOTCAR( ADRESS(4, 4) )
508 ! NOMCAS = MOTCAR( ADRESS(4, 5) )
509  wac_files(waccli)%NAME=motcar( adress(4,6) )
510  IF(wac_files(waccli)%NAME(1:1).EQ.' ') THEN
511  WRITE(lu,*) 'THE FOLLOWING KEYWORD IS MANDATORY:'
512  WRITE(lu,*) 'BOUNDARY CONDITIONS FILE '//
513  & '(FICHIER DES CONDITIONS AUX LIMITES)'
514  CALL plante(1)
515  stop
516  ENDIF
517  wac_files(wacfon)%NAME=motcar( adress(4,7) )
518  wac_files(wacres)%NAME=motcar( adress(4,8) )
519  wac_files(wacleo)%NAME=motcar( adress(4,9) )
520  wac_files(wacpre)%NAME=motcar( adress(4,10) )
521  wac_files(wacrbi)%NAME=motcar( adress(4,11) )
522  wac_files(waccob)%NAME=motcar( adress(4,12) )
523  wac_files(waccof)%NAME=motcar( adress(4,13) )
524  wac_files(wacbi1)%NAME=motcar( adress(4,14) )
525  wac_files(wacfo1)%NAME=motcar( adress(4,15) )
526  wac_files(leowxy)%NAME=motcar( adress(4,23) )
527  wac_files(leoixy)%NAME=motcar( adress(4,24) )
528  wac_files(impspe)%NAME=motcar( adress(4,25) )
529  IF( wac_files(wacrbi)%NAME.NE.' ') THEN
530 ! ONE WANTS TO HAVE A GLOBAL RESULT
531  glob=.true.
532  IF (inclus(coupling,'TOMAWAC')) THEN
533  IF(abs(nit*dt-cpl_wac_data%NIT_TEL*cpl_wac_data%DT_TEL)
534  & .GE.1d-6)THEN
535  WRITE(lu,*) nit,dt,cpl_wac_data%NIT_TEL,cpl_wac_data%DT_TEL
536  WRITE(lu,*)'WHEN ONE WANTS A GLOBAL RESULT WITH COUPLING'
537  WRITE(lu,*)'DURATION OF TELEMAC AND TOMAWAC MUST',
538  & 'BE THE SAME'
539  CALL plante(1)
540  ENDIF
541  IF(mod(nit,cpl_wac_data%PERCOU_WAC).NE.0)THEN
542  WRITE(lu,*)'WHEN ONE WANTS A GLOBAL RESULT WITH COUPLING'
543  WRITE(lu,*)'NUMBER OF TIME STEP OF TOMAWAC MUST',
544  & 'BE A MULTIPLE OF COUPLING PERIOD'
545  CALL plante(1)
546  ENDIF
547  ENDIF
548  ELSE
549  glob=.false.
550  ENDIF
551  IF( wac_files(wacpre)%NAME.NE.' ') THEN
552 ! ONE WANTS TO HAVE A GLOBAL RESULT
553  suit=.true.
554  ELSE
555  suit=.false.
556  ENDIF
557 !
558 ! FROM 26 TO 28 : FOR CRAY, NOT USEFUL HERE
559 !
560  wac_files(wacveb)%NAME=motcar( adress(4,31) )
561  wac_files(wacvef)%NAME=motcar( adress(4,32) )
562  wac_files(wacref)%NAME=motcar( adress(4,34) )
563  wac_files(wacmab)%NAME=motcar( adress(4,35) )
564  wac_files(wacmaf)%NAME=motcar( adress(4,36) )
565  equa = 'TOMAWAC-COWADIS'
566 !BD_INCKA FILE FORMATS
567 ! GEOMETRY FILE
568  wac_files(wacgeo)%FMT = motcar( adress(4,39) )(1:8)
569  CALL majus(wac_files(wacgeo)%FMT)
570 ! RESULTS FILE FORMAT
571  wac_files(wacres)%FMT = motcar( adress(4,17) )(1:8)
572  CALL majus(wac_files(wacres)%FMT)
573 ! INITIAL RESULTS FILE FORMAT (< PREVIOUS COMPUTATION)
574 ! SEDIMENT...
575  wac_files(wacpre)%FMT = motcar( adress(4,41) )(1:8)
576  CALL majus(wac_files(wacpre)%FMT)
577 ! REFERENCE FILE FORMAT
578  wac_files(wacref)%FMT = motcar( adress(4,42) )(1:8)
579  CALL majus(wac_files(wacref)%FMT)
580 ! BINARY FILE 1 FORMAT
581  wac_files(wacbi1)%FMT = motcar( adress(4,43) )(1:8)
582  CALL majus(wac_files(wacbi1)%FMT)
583 ! SPECTRAL FILE FORMAT
584  wac_files(wacleo)%FMT = motcar( adress(4,44) )(1:8)
585  CALL majus(wac_files(wacleo)%FMT)
586 ! GLOBAL RESULT FILE FORMAT
587  wac_files(wacrbi)%FMT = motcar( adress(4,51) )(1:8)
588  CALL majus(wac_files(wacrbi)%FMT)
589 ! BINARY CURRENTS FILE FORMAT
590  wac_files(waccob)%FMT = motcar( adress(4,40) )(1:8)
591  CALL majus(wac_files(waccob)%FMT)
592 ! BINARY WINDS FILE FORMAT
593  wac_files(wacveb)%FMT = motcar( adress(4,46) )(1:8)
594  CALL majus(wac_files(wacveb)%FMT)
595 ! BINARY TIDAL WATER FILE FORMAT
596  wac_files(wacmab)%FMT = motcar( adress(4,47) )(1:8)
597  CALL majus(wac_files(wacmab)%FMT)
598 ! IMPOSED SPECTRA FILE FORMAT
599  wac_files(impspe)%FMT = motcar( adress(4,26) )(1:8)
600  CALL majus(wac_files(impspe)%FMT)
601 
602 !
603 ! NAMES OF VARIABLES
604 !
605  nameu =motcar( adress(4,45) )(1:32)
606  namev =motcar( adress(4,45)+1 )(1:32)
607  namewx=motcar( adress(4,45)+2 )(1:32)
608  namewy=motcar( adress(4,45)+3 )(1:32)
609  nameh =motcar( adress(4,45)+4 )(1:32)
610 !
611  wac_files(wacspe)%NAME=motcar( adress(4,50) )
612 !
613 ! CORRECTS OR COMPUTES OTHER PARAMETERS FROM THOSE THAT
614 ! HAVE JUST BEEN READ
615 !
616  IF(cousta.OR.maree) THEN
617  couran=.true.
618  ELSE
619  couran=.false.
620  ENDIF
621  IF(.NOT.vent.AND.svent.NE.0) THEN
622  WRITE(lu,*)
623  & 'INCOMPATIBILITY OF KEY-WORDS CONCERNING WIND => NO
624  & WIND'
625  svent=0
626  ENDIF
627  IF(trigo) THEN
628  teta1 = pisur2-teta1
629  teta2 = pisur2-teta2
634  ENDIF
635  IF(cimpli.LT.0.OR.cimpli.GT.1) THEN
636  WRITE(lu,*) 'INCOMPATIBILITY OF IMPLICITATION COEFFICIENT'
637  WRITE(lu,*) 'VALUE READ = ',cimpli
638  WRITE(lu,*) 'WE TAKE THE DEFAULT VALUE CIMPLI=0.5'
639  cimpli=0.5d0
640  ENDIF
641 ! NEW SOURCE TERMS
642  IF(.NOT.proinf.AND.strif.EQ.3) THEN
643  WRITE(lu,*) 'INCOMPATIBILITY OF DEPTH AND'
644  WRITE(lu,*) 'NON-LINEAR TRANSFER TERM'
645  CALL plante(1)
646  stop
647  ENDIF
648  IF(cbaj.EQ.1) THEN
649  limit=3
650 ! IT used to include non linear transfert =1
651 ! Now user has to set a non linear Law.
652 ! STRIF=1
653  svent=1
654  smout=1
655  cmout1=2.1d0
656  cmout2=0.4d0
657  alpha=0.0095d0
658  ENDIF
659 !
660  IF ((npleo.EQ.0.AND.namwxy(1:1).EQ.' ' ).AND.
661  & (namleo(1:1).NE.' '.OR.namspe(1:1).NE.' ')) THEN
662 
663  WRITE(lu,*) ''
664  WRITE(lu,*) '********************************************'
665  WRITE(lu,*) ' A PUNCTUAL NAME FILE HAS BEEN DEFINED '
666  WRITE(lu,*) ' BUT NO COORDINATE HAS BEEN DEFINED '
667  WRITE(lu,*) '********************************************'
668  CALL plante(1)
669  ENDIF
670  IF ((npleo.GT.0.OR.namwxy(1:1).NE.' ' ) .AND.
671  & namleo(1:1).EQ.' '.AND.namspe(1:1).EQ.' ') THEN
672  WRITE(lu,*) ''
673  WRITE(lu,*) '********************************************'
674  WRITE(lu,*) ' SOME COORDINATES HAVE BEEN DEFINED '
675  WRITE(lu,*) ' BUT NO PUNCTUAL NAME FILE HAS BEEN DEFINED'
676  WRITE(lu,*) '********************************************'
677  CALL plante(1)
678  ENDIF
679 
680 !
681 !
682 !-----------------------------------------------------------------------
683 ! NAME OF THE VARIABLES FOR THE RESULTS AND GEOMETRY FILES:
684 !-----------------------------------------------------------------------
685 !
686 ! LOGICAL ARRAY FOR OUTPUT
687 !
688  CALL nomvar_tomawac(texte,mnemo,maxvar)
689 !
690 !$DC$ BUG : ARRAYS MNEMO AND SORLEO OF SIZE MAXVAR
691 ! MUCH LESS THAN 100 !
692  CALL sortie(sort2d , mnemo , maxvar , sorleo )
693 !
694 !.....IF NO WIND, THERE SHOULD BE NO INFORMATION WRITTEN ABOUT WINDS
695  IF (.NOT.vent) THEN
696  sorleo( 9)=.false.
697  sorleo(10)=.false.
698  sorleo(25)=.false.
699  ENDIF
700 !.....IF NO CURRENT, THERE SHOULD BE NO INFORMATION WRITTEN ABOUT COURANT
701  IF (.NOT.couran.AND..NOT.inclus(coupling,'TOMAWAC')) THEN
702  sorleo(7)=.false.
703  sorleo(8)=.false.
704  IF (sdscu.NE.0) THEN
705 ! AND NO DISSIPATION COEFFICIENT FOR STRONG CURRENT
706  WRITE(lu,*) '*****************************************'
707  WRITE(lu,*) ' NO DISSIPATION FOR STRONG CURRENT'
708  WRITE(lu,*) '*****************************************'
709  sdscu = 0
710  ENDIF
711  ENDIF
712 !
713 !.....IF INFINITE DEPTH, THE RADIATION STRESSES ARE NOT COMPUTED
714  IF (proinf) THEN
715  IF (sorleo(11) .OR. sorleo(12) .OR. sorleo(13) .OR.
716  & sorleo(14) .OR. sorleo(15) ) THEN
717  WRITE(lu,*) '*****************************************'
718  WRITE(lu,*) ' RADIATION STRESSES ARE NOT COMPUTED '
719  WRITE(lu,*) ' OVER INFINITE WATER DEPTHS '
720  WRITE(lu,*) '******************************************'
721  DO k=11,15
722  sorleo(k) = .false.
723  ENDDO
724  ENDIF
725 ! LA HAUTEUR D'EAU ET LA PROFONDEUR NE DOIVENT PAS ETRE CALCULEE
726 ! WATER DEPTH AND DEPTH MUST NOT BE COMPUTED
727  sorleo(5) = .false.
728  sorleo(6) = .false.
729  ENDIF
730 !
731  DO k=1,maxvar
732  sorimp(k)=.false.
733  ENDDO
734 !
735 !-----------------------------------------------------------------------
736 !
737  RETURN
738  END
subroutine tomawac_constants
character(len=path_len), pointer namwxy
character(len=8), pointer fmtgeo
logical, dimension(maxvar) sorimp
character(len=8), pointer fmtpre
character(len=8), pointer fmtveb
character(len=path_len), pointer namixy
double precision, dimension(:), allocatable yleo
integer, parameter maxvar
character(len=8), pointer fmtvef
character(len=80) titcas
subroutine nomvar_tomawac(TEXTE, MNEMO, MAXVAR)
Definition: nomvar_tomawac.f:7
subroutine read_submit(FILES, NFILES, SUBMIT, NMOT)
Definition: read_submit.f:7
integer, parameter maxlu_wac
character(len=8), pointer fmtmab
character(len=8), pointer fmtres
character(len=32), dimension(maxvar) texte
integer, parameter maxkeyword
character(len=path_len), pointer nammaf
character(len=path_len), pointer namspe
character(len=8), pointer fmtref
type(cpl_wac_data_obj) cpl_wac_data
character(len=path_len), pointer namvef
character(len=8), pointer fmtbi1
character(len=path_len), pointer namcob
logical function inclus(C1, C2)
Definition: inclus.f:7
subroutine sortie(CHAINE, MNEMO, NBRE, SORLEO)
Definition: sortie.f:7
subroutine damocle(ADRESS, DIMENS, NMAX, DOC, LLNG, LLU, MOTINT, MOTREA, MOTLOG, MOTCAR, MOTCLE, TROUVE, NFICMO, NFICDA, GESTD, MOTATT)
Definition: damocle.f:9
logical, dimension(maxvar) sorleo
character(len=8), pointer fmtcob
character(len=path_len), pointer namfon
character(len=8), pointer fmtmaf
character(len=8), pointer fmtrbi
character(len=8), pointer fmtleo
character(len=8), pointer fmtcof
character(len=path_len), pointer namcof
integer, dimension(:), allocatable noleo
subroutine lecdon_tomawac(FILE_DESC, PATH, NCAR, CAS_FILE, DICO_FILE)
Definition: lecdon_tomawac.f:7
character(len=72) sort2d
double precision, dimension(:), allocatable xleo
character(len=path_len), pointer namveb
character(len=path_len), pointer nammab
double precision, dimension(:), pointer x
subroutine majus(CHAINE)
Definition: majus.f:7
character(len=path_len), pointer namres
character(len=path_len), target coupling
character(len=path_len), pointer namleo
double precision, target dt
type(bief_file), dimension(maxlu_wac), target wac_files
Definition: bief.f:3