IMEX_SfloW2D  0.9
Shallowwatergranularflowmodel
inpout_2d.f90
Go to the documentation of this file.
1 !********************************************************************************
3 !
6 !
10 !
11 !********************************************************************************
12 
13 MODULE inpout_2d
14 
15  ! -- Variables for the namelist RUN_PARAMETERS
16  USE parameters_2d, ONLY : t_start , t_end , dt_output , &
17  topography_function_flag, topography_demfile
18 
19  USE solver_2d, ONLY : verbose_level
20 
21  ! -- Variables for the namelist NEWRUN_PARAMETERS
22  USE geometry_2d, ONLY : x0 , y0 , comp_cells_x , comp_cells_y , cell_size
23  USE geometry_2d, ONLY : topography_profile , n_topography_profile_x , &
24  n_topography_profile_y
25  USE init_2d, ONLY : riemann_interface
26  USE parameters_2d, ONLY : temperature_flag , solid_transport_flag , &
27  riemann_flag , rheology_flag
28  USE parameters_2d, ONLY : source_flag
29 
30  ! -- Variables for the namelist INITIAL_CONDITIONS
31  USE parameters_2d, ONLY : released_volume , x_release , y_release
32  USE parameters_2d, ONLY : velocity_mod_release , velocity_ang_release
33  USE parameters_2d, ONLY : xs_init , xs_ambient , sed_vol_perc
34  USE parameters_2d, ONLY : t_init , t_ambient
35 
36  ! -- Variables for the namelist LEFT_STATE
37  USE init_2d, ONLY : hb_w , u_w , v_w , xs_w , t_w
38 
39  ! -- Variables for the namelist RIGHT_STATE
40  USE init_2d, ONLY : hb_e , u_e , v_e , xs_e , t_e
41 
42  ! -- Variables for the namelists LEFT/RIGHT_BOUNDARY_CONDITIONS
43  USE parameters_2d, ONLY : bc
44 
45  ! -- Variables for the namelist NUMERIC_PARAMETERS
46  USE parameters_2d, ONLY : solver_scheme, dt0 , max_dt , cfl, limiter , theta, &
47  reconstr_variables , reconstr_coeff , interfaces_relaxation , n_rk
48 
49  ! -- Variables for the namelist EXPL_TERMS_PARAMETERS
50  USE constitutive_2d, ONLY : grav
51  USE parameters_2d, ONLY : x_source , y_source , r_source , vfr_source , &
52  t_source
53 
54  ! -- Variables for the namelist TEMPERATURE_PARAMETERS
55  USE constitutive_2d, ONLY : rho , emissivity , exp_area_fract , enne , emme , &
56  atm_heat_transf_coeff , thermal_conductivity , t_env , t_ground , c_p
57 
58  ! -- Variables for the namelist RHEOLOGY_PARAMETERS
59  USE parameters_2d, ONLY : rheology_model
60  USE constitutive_2d, ONLY : mu , xi , tau , nu_ref , visc_par , t_ref
61  USE constitutive_2d, ONLY : alpha2 , beta2 , alpha1 , beta1 , kappa , n_td , &
62  gamma_w , gamma_s
63 
64  IMPLICIT NONE
65 
66  CHARACTER(LEN=40) :: run_name
67  CHARACTER(LEN=40) :: bak_name
68  CHARACTER(LEN=40) :: input_file
69  CHARACTER(LEN=40) :: output_file
70  CHARACTER(LEN=40) :: restart_file
71  CHARACTER(LEN=40) :: probes_file
72  CHARACTER(LEN=40) :: output_file_2d
73  CHARACTER(LEN=40) :: output_esri_file
74  CHARACTER(LEN=40) :: runout_file
75 
76  INTEGER, PARAMETER :: input_unit = 7
77  INTEGER, PARAMETER :: backup_unit = 8
78  INTEGER, PARAMETER :: output_unit = 9
79  INTEGER, PARAMETER :: restart_unit = 10
80  INTEGER, PARAMETER :: probes_unit = 11
81  INTEGER, PARAMETER :: output_unit_2d = 12
82  INTEGER, PARAMETER :: output_esri_unit = 13
83  INTEGER, PARAMETER :: dem_esri_unit = 14
84  INTEGER, PARAMETER :: runout_unit = 15
85  INTEGER, PARAMETER :: dakota_unit = 16
86 
88  INTEGER :: output_idx
89 
94  LOGICAL :: restart
95 
100  LOGICAL :: output_esri_flag
101 
106  LOGICAL :: output_phys_flag
107 
112  LOGICAL :: output_cons_flag
113 
118  LOGICAL :: output_runout_flag
119 
120  ! -- Variables for the namelists WEST_BOUNDARY_CONDITIONS
121  TYPE(bc) :: hb_bcw , u_bcw , v_bcw , hu_bcw , hv_bcw , xs_bcw , t_bcw
122 
123  ! -- Variables for the namelists EAST_BOUNDARY_CONDITIONS
124  TYPE(bc) :: hb_bce , u_bce , v_bce , hu_bce , hv_bce , xs_bce , t_bce
125 
126  ! -- Variables for the namelists SOUTH_BOUNDARY_CONDITIONS
127  TYPE(bc) :: hb_bcs , u_bcs , v_bcs , hu_bcs , hv_bcs , xs_bcs , t_bcs
128 
129  ! -- Variables for the namelists NORTH_BOUNDARY_CONDITIONS
130  TYPE(bc) :: hb_bcn , u_bcn , v_bcn , hu_bcn , hv_bcn , xs_bcn , t_bcn
131 
132 
133  ! parameters to read a dem file
134  INTEGER :: ncols, nrows, nodata_value
135 
136  REAL*8 :: xllcorner, yllcorner, cellsize
137 
138  LOGICAL :: write_first_q
139 
140  INTEGER :: n_probes
141 
142  REAL*8, ALLOCATABLE :: probes_coords(:,:)
143 
144  REAL*8 :: dt_runout
145 
146  REAL*8, ALLOCATABLE :: h_old(:,:)
147 
148  REAL*8 :: x0_runout, y0_runout , init_runout , runout_mom_stop
149 
150  namelist / run_parameters / run_name , restart , topography_demfile , &
151  t_start , t_end , dt_output , output_cons_flag , output_esri_flag , &
152  output_phys_flag , output_runout_flag , verbose_level
153 
154  namelist / restart_parameters / restart_file , xs_init , xs_ambient , t_init ,&
155  t_ambient , sed_vol_perc
156 
157  namelist / newrun_parameters / x0 , y0 , comp_cells_x , comp_cells_y , &
158  cell_size , temperature_flag , source_flag , solid_transport_flag , &
159  rheology_flag , riemann_flag
160 
161  namelist / initial_conditions / released_volume , x_release , y_release , &
162  velocity_mod_release , velocity_ang_release , t_init , t_ambient
163 
164  namelist / left_state / riemann_interface , hb_w , u_w , v_w , xs_w , t_w
165 
166  namelist / right_state / hb_e , u_e , v_e , xs_e , t_e
167 
168  namelist / west_boundary_conditions / hb_bcw , u_bcw , v_bcw , hu_bcw , &
169  hv_bcw , xs_bcw , t_bcw
170 
171  namelist / east_boundary_conditions / hb_bce , u_bce , v_bce , hu_bce , &
172  hv_bce , xs_bce , t_bce
173 
174  namelist / south_boundary_conditions / hb_bcs , u_bcs , v_bcs , hu_bcs , &
175  hv_bcs , xs_bcs , t_bcs
176 
177  namelist / north_boundary_conditions / hb_bcn , u_bcn , v_bcn , hu_bcn , &
178  hv_bcn , xs_bcn , t_bcn
179 
180  namelist / numeric_parameters / solver_scheme, dt0 , max_dt , cfl, limiter , &
181  theta , reconstr_variables , reconstr_coeff , interfaces_relaxation , n_rk
182 
183  namelist / expl_terms_parameters / grav , x_source , y_source , r_source , &
184  vfr_source , t_source
185 
186  namelist / temperature_parameters / emissivity , atm_heat_transf_coeff , &
187  thermal_conductivity , exp_area_fract , c_p , enne , emme , t_env , &
188  t_ground , rho
189 
190  namelist / rheology_parameters / rheology_model , mu , xi , tau , nu_ref , &
191  visc_par , t_ref , alpha2 , beta2 , alpha1 , beta1 , kappa , n_td , &
192  gamma_w , gamma_s
193 
194  namelist / runout_parameters / x0_runout , y0_runout , dt_runout , &
195  runout_mom_stop
196 
197 CONTAINS
198 
199  !******************************************************************************
201  !
205  !
209  !
210  !******************************************************************************
211 
212  SUBROUTINE init_param
213 
214  USE parameters_2d , ONLY : n_vars
215 
216  IMPLICIT none
217 
218  LOGICAL :: lexist
219 
220  INTEGER :: j , k
221 
222  !-- Inizialization of the Variables for the namelist RUN_PARAMETERS
223  run_name = 'default'
224  restart = .false.
225  topography_function_flag=.false.
226  topography_demfile=.false.
227  t_start = 0.0
228  t_end = 5.0d-4
229  dt_output = 1.d-4
230  output_cons_flag = .true.
231  output_esri_flag = .true.
232  output_phys_flag = .true.
233  output_runout_flag = .false.
234  verbose_level = 0
235 
236  !-- Inizialization of the Variables for the namelist restart parameters
237  restart_file = ''
238  xs_init = 0.d0
239  xs_ambient = 0.d0
240  t_init = 0.d0
241  t_ambient = 0.d0
242  sed_vol_perc = 0.d0
243 
244  !-- Inizialization of the Variables for the namelist newrun_parameters
245  x0 = 0.d0
246  comp_cells_x = 400
247  y0 = 0.d0
248  comp_cells_y = 1
249  cell_size = 2.5d-3
250  temperature_flag = .false.
251  source_flag = .false.
252  rheology_flag = .false.
253  riemann_flag =.true.
254  solid_transport_flag = .false.
255 
256  !-- Inizialization of the Variables for the namelist left_state
257  riemann_interface = 0.5d0
258  hb_w = 2.222d0
259  u_w = 0.8246d0
260  v_w = 0.d0
261  xs_w = 0.5d0
262  t_w = -1.d0
263 
264  !-- Inizialization of the Variables for the namelist right_state
265  hb_e = -1.0d0
266  u_e = -1.6359d0
267  u_e = 0.d0
268  xs_e = 0.d0
269  t_e = -1.d0
270 
271  !-- Inizialization of the Variables for the namelist west boundary conditions
272  hb_bcw%flag = 1
273  hb_bcw%value = 0.d0
274 
275  u_bcw%flag = 1
276  u_bcw%value = 0.d0
277 
278  v_bcw%flag = 1
279  v_bcw%value = 0.d0
280 
281  xs_bcw%flag = 1
282  xs_bcw%value = 0.d0
283 
284  !-- Inizialization of the Variables for the namelist east boundary conditions
285  hb_bce%flag = 1
286  hb_bce%value = 0.d0
287 
288  u_bce%flag = 1
289  u_bce%value = 0.d0
290 
291  v_bce%flag = 1
292  v_bce%value = 0.d0
293 
294  xs_bce%flag = 1
295  xs_bce%value = 0.d0
296 
297  !-- Inizialization of the Variables for the namelist south boundary conditions
298  hb_bcs%flag = 1
299  hb_bcs%value = 0.d0
300 
301  u_bcs%flag = 1
302  u_bcs%value = 0.d0
303 
304  v_bcs%flag = 1
305  v_bcs%value = 0.d0
306 
307  xs_bcs%flag = 1
308  xs_bcs%value = 0.d0
309 
310  !-- Inizialization of the Variables for the namelist north boundary conditions
311  hb_bcn%flag = 1
312  hb_bcn%value = 0.d0
313 
314  u_bcn%flag = 1
315  u_bcn%value = 0.d0
316 
317  v_bcn%flag = 1
318  v_bcn%value = 0.d0
319 
320  xs_bcn%flag = 1
321  xs_bcn%value = 0.d0
322 
323  !-- Inizialization of the Variables for the namelist NUMERIC_PARAMETERS
324  dt0 = 1.d-4
325  max_dt = 1.d-3
326  solver_scheme = 'KT'
327  n_rk = 1
328  cfl = 0.24
329  limiter(1:n_vars) = 0
330  theta=1.0
331  reconstr_variables = 'phys'
332  reconstr_coeff = 1.0
333 
334  !-- Inizialization of the Variables for the namelist EXPL_TERMS_PARAMETERS
335  grav = 9.81d0
336  x_source = 0.d0
337  y_source = 0.d0
338  r_source = 0.d0
339  vfr_source = 0.d0
340  t_source = 0.d0
341 
342  !-- Inizialization of the Variables for the namelist TEMPERATURE_PARAMETERS
343  exp_area_fract = 0.5d0
344  emissivity = 0.0d0 ! no radiation to atmosphere
345  atm_heat_transf_coeff = 0.0d0 ! no convection to atmosphere
346  thermal_conductivity = 0.0d0 ! no conduction to ground
347  enne = 4.0d0
348  emme = 12.d0
349  t_env = 300.0d0
350  t_ground = 1200.0d0
351  c_p = 1200.d0
352  rho = 0.d0
353 
354  !-- Inizialization of the Variables for the namelist RHEOLOGY_PARAMETERS
355  rheology_model = 0
356  nu_ref = 0.0d0
357  mu = 0.d0
358  xi = 0.d0
359  tau = 0.d0
360  t_ref = 0.d0
361  visc_par = 0.0d0
362 
363  !-- Inizialization of the Variables for the namelist RUNOUT_PARAMETERS
364  x0_runout = -1
365  y0_runout = -1
366  dt_runout = 60
367  runout_mom_stop = 0.d0
368 
369  !-------------- Check if input file exists ----------------------------------
370  input_file = 'IMEX_SfloW2D.inp'
371 
372  INQUIRE (file=input_file,exist=lexist)
373 
374  IF (lexist .EQV. .false.) THEN
375 
376  OPEN(input_unit,file=input_file,status='NEW')
377 
378  WRITE(input_unit, run_parameters )
379  WRITE(input_unit, newrun_parameters )
380  WRITE(input_unit, left_state )
381  WRITE(input_unit, right_state )
382  WRITE(input_unit, numeric_parameters )
383  WRITE(input_unit, west_boundary_conditions )
384  WRITE(input_unit, east_boundary_conditions )
385  WRITE(input_unit, south_boundary_conditions )
386  WRITE(input_unit, north_boundary_conditions )
387  WRITE(input_unit, expl_terms_parameters )
388 
389  n_topography_profile_x = 2
390  n_topography_profile_y = 2
391 
392  ALLOCATE( topography_profile( 2 , n_topography_profile_x , &
393  n_topography_profile_y) )
394 
395  topography_profile(1,1,1) = 0.d0
396  topography_profile(2,1,1) = 0.d0
397  topography_profile(3,1,1) = 1.d0
398 
399  topography_profile(1,1,2) = 0.d0
400  topography_profile(2,1,2) = 1.d0
401  topography_profile(3,1,2) = 1.d0
402 
403  topography_profile(1,2,1) = 1.d0
404  topography_profile(2,2,1) = 0.d0
405  topography_profile(3,2,1) = 0.5d0
406 
407  topography_profile(1,2,2) = 1.d0
408  topography_profile(2,2,2) = 1.d0
409  topography_profile(3,2,2) = 0.5d0
410 
411 
412 
413  WRITE(input_unit,*) '''TOPOGRAPHY_PROFILE'''
414  WRITE(input_unit,*) n_topography_profile_x
415  WRITE(input_unit,*) n_topography_profile_y
416 
417  DO j = 1, n_topography_profile_x
418 
419  DO k = 1, n_topography_profile_y
420 
421  WRITE(input_unit,108) topography_profile(1:3,j,k)
422 
423 108 FORMAT(3(1x,e14.7))
424 
425  ENDDO
426 
427  END DO
428 
429  CLOSE(input_unit)
430 
431  WRITE(*,*) 'Input file IMEX_SfloW2D.inp not found'
432  WRITE(*,*) 'A new one with default values has been created'
433  stop
434 
435  ELSE
436 
437  END IF
438 
439  ! output file index
440  output_idx = 0
441 
442 
443  ! -------------- Initialize values for checks during input reading ----------
444  hb_bcw%flag = -1
445  u_bcw%flag = -1
446  v_bcw%flag = -1
447  hu_bcw%flag = -1
448  hv_bcw%flag = -1
449  xs_bcw%flag = -1
450  t_bcw%flag = -1
451 
452  hb_bce%flag = -1
453  u_bce%flag = -1
454  v_bce%flag = -1
455  hu_bce%flag = -1
456  hv_bce%flag = -1
457  xs_bce%flag = -1
458  t_bce%flag = -1
459 
460  hb_bcs%flag = -1
461  u_bcs%flag = -1
462  v_bcs%flag = -1
463  hu_bcs%flag = -1
464  hv_bcs%flag = -1
465  xs_bcs%flag = -1
466  t_bcs%flag = -1
467 
468  hb_bcn%flag = -1
469  u_bcn%flag = -1
470  v_bcn%flag = -1
471  hu_bcn%flag = -1
472  hv_bcn%flag = -1
473  xs_bcn%flag = -1
474  t_bcn%flag = -1
475 
476  sed_vol_perc = -1.d0
477 
478  rheology_model = -1
479  mu = -1
480  xi = -1
481  tau = -1
482  nu_ref = -1
483  visc_par = -1
484  t_ref = -1
485 
486  alpha2 = -1
487  beta2 = -1
488  alpha1 = -1
489  beta1 = -1
490  kappa = -1
491  n_td = -1
492  gamma_w = -1
493  gamma_s = -1
494 
495  exp_area_fract = -1.d0
496  emissivity = -1.d0
497  atm_heat_transf_coeff = -1.d0
498  thermal_conductivity = -1.d0
499  enne = -1.d0
500  emme = -1.d0
501  t_env = -1.d0
502  t_ground = -1.d0
503  c_p = -1.d0
504  rho = -1.d0
505 
506  grav = -1.d0
507 
508  x0_runout = -1.d0
509  y0_runout = -1.d0
510 
511 
512  END SUBROUTINE init_param
513 
514  !******************************************************************************
516  !
521  !
525  !
526  !******************************************************************************
527 
528  SUBROUTINE read_param
529 
530  USE parameters_2d, ONLY : n_vars , n_eqns
531  USE parameters_2d, ONLY : limiter
532  USE parameters_2d, ONLY : bcw , bce , bcs , bcn
533 
534  IMPLICIT none
535 
536  REAL*8 :: max_cfl
537 
538  LOGICAL :: tend1
539  CHARACTER(LEN=80) :: card
540 
541  INTEGER :: j,k
542 
543  INTEGER :: dot_idx
544 
545  CHARACTER(LEN=3) :: check_file
546 
547  LOGICAL :: lexist
548 
549  CHARACTER(LEN=15) :: chara
550 
551  INTEGER :: ios
552 
553  OPEN(input_unit,file=input_file,status='old')
554 
555 
556  ! ------- READ run_parameters NAMELIST -----------------------------------
557  READ(input_unit, run_parameters,iostat=ios )
558 
559  IF ( ios .NE. 0 ) THEN
560 
561  WRITE(*,*) 'IOSTAT=',ios
562  WRITE(*,*) 'ERROR: problem with namelist RUN_PARAMETERS'
563  WRITE(*,*) 'Please check the input file'
564  stop
565 
566  ELSE
567 
568  rewind(input_unit)
569 
570  END IF
571 
572  ! ------- READ newrun_parameters NAMELIST --------------------------------
573  READ(input_unit,newrun_parameters,iostat=ios)
574 
575  IF ( ios .NE. 0 ) THEN
576 
577  WRITE(*,*) 'IOSTAT=',ios
578  WRITE(*,*) 'ERROR: problem with namelist NEWRUN_PARAMETERS'
579  WRITE(*,*) 'Please check the input file'
580  stop
581 
582  ELSE
583 
584  rewind(input_unit)
585 
586  END IF
587 
588  IF ( ( comp_cells_x .EQ. 1 ) .OR. ( comp_cells_y .EQ. 1 ) ) THEN
589 
590  WRITE(*,*) '----- 1D SIMULATION -----'
591 
592  ELSE
593 
594  WRITE(*,*) '----- 2D SIMULATION -----'
595 
596  END IF
597 
598  IF ( solid_transport_flag ) THEN
599 
600  IF ( temperature_flag ) THEN
601 
602  n_vars = 5
603  n_eqns = 5
604 
605  ELSE
606 
607  n_vars = 4
608  n_eqns = 4
609 
610  END IF
611 
612  ELSE
613 
614  IF ( temperature_flag ) THEN
615 
616  n_vars = 4
617  n_eqns = 4
618 
619  ELSE
620 
621  n_vars = 3
622  n_eqns = 3
623 
624  END IF
625 
626  END IF
627 
628  ALLOCATE( bcw(n_vars) , bce(n_vars) , bcs(n_vars) , bcn(n_vars) )
629 
630  IF ( restart ) THEN
631 
632  READ(input_unit,restart_parameters,iostat=ios)
633 
634  IF ( ios .NE. 0 ) THEN
635 
636  WRITE(*,*) 'IOSTAT=',ios
637  WRITE(*,*) 'ERROR: problem with namelist RESTART_PARAMETERS'
638  WRITE(*,*) 'Please check the input file'
639  stop
640 
641  ELSE
642 
643  IF ( ( sed_vol_perc .LT. -1.d0 ) .OR. ( sed_vol_perc .GT. 100.d0 ) ) &
644  THEN
645 
646  WRITE(*,*) 'ERROR: problem with namelist RESTART_PARAMETERS'
647  WRITE(*,*) 'SED_VOL_PERC =' , sed_vol_perc
648  stop
649 
650  END IF
651 
652  xs_init = 1.d-2 * sed_vol_perc
653  xs_ambient = 1.d-2 * sed_vol_perc
654 
655  rewind(input_unit)
656 
657  END IF
658 
659  IF ( temperature_flag .AND. restart ) THEN
660 
661  dot_idx = scan(restart_file, ".", .true.)
662 
663  check_file = restart_file(dot_idx+1:dot_idx+3)
664 
665  IF ( ( check_file .EQ. 'asc' ) .AND. &
666  ( t_init*t_ambient .EQ. 0.d0 ) ) THEN
667 
668  WRITE(*,*) 'T_init=',t_init
669  WRITE(*,*) 'T_ambient=',t_ambient
670  WRITE(*,*) 'Add the variables to the namelist RESTART_PARAMETERS'
671  stop
672 
673  END IF
674 
675  END IF
676 
677  ELSE
678 
679  IF ( riemann_flag ) THEN
680 
681  READ(input_unit,left_state,iostat=ios)
682 
683  IF ( ios .NE. 0 ) THEN
684 
685  WRITE(*,*) 'IOSTAT=',ios
686  WRITE(*,*) 'ERROR: problem with namelist LEFT_PARAMETERS'
687  WRITE(*,*) 'Please check the input file'
688  stop
689 
690  ELSE
691 
692  rewind(input_unit)
693 
694  IF ( ( temperature_flag ) .AND. ( t_w .EQ. -1.d0 ) ) THEN
695 
696  WRITE(*,*) 'ERROR: problem with namelist LEFT_PARAMETERS'
697  WRITE(*,*) 'Initial temperature T_W not defined'
698  stop
699 
700  END IF
701 
702  END IF
703 
704  READ(input_unit,right_state,iostat=ios)
705 
706  IF ( ios .NE. 0 ) THEN
707 
708  WRITE(*,*) 'IOSTAT=',ios
709  WRITE(*,*) 'ERROR: problem with namelist RIGHT_PARAMETERS'
710  WRITE(*,*) 'Please check the input file'
711  stop
712 
713  ELSE
714 
715  rewind(input_unit)
716 
717  IF ( ( temperature_flag ) .AND. ( t_e .EQ. -1.d0 ) ) THEN
718 
719  WRITE(*,*) 'ERROR: problem with namelist RIGHT_PARAMETERS'
720  WRITE(*,*) 'Initial temperature T_E not defined'
721  stop
722 
723  END IF
724 
725 
726  END IF
727 
728  ELSE
729 
730  READ(input_unit,initial_conditions,iostat=ios)
731 
732  IF ( ios .NE. 0 ) THEN
733 
734  WRITE(*,*) 'IOSTAT=',ios
735  WRITE(*,*) 'ERROR: problem with namelist INITIAL_CONDITIONS'
736  WRITE(*,*) 'Please check the input file'
737  stop
738 
739  ELSE
740 
741  rewind(input_unit)
742 
743  END IF
744 
745  IF ( ( temperature_flag ) .AND. ( t_init*t_ambient .EQ. 0.d0 ) ) THEN
746 
747  WRITE(*,*) 'T_init=',t_init
748  WRITE(*,*) 'T_ambient=',t_ambient
749  WRITE(*,*) 'Add the two variables to namelist INITIAL_CONDITIONS'
750  stop
751 
752  END IF
753 
754  END IF
755 
756  END IF
757 
758  ! ------- READ numeric_parameters NAMELIST ----------------------------------
759 
760  READ(input_unit,numeric_parameters)
761 
762  IF ( ios .NE. 0 ) THEN
763 
764  WRITE(*,*) 'IOSTAT=',ios
765  WRITE(*,*) 'ERROR: problem with namelist NUMERIC_PARAMETERS'
766  WRITE(*,*) 'Please check the input file'
767  stop
768 
769  ELSE
770 
771  rewind(input_unit)
772 
773  END IF
774 
775  IF ( ( solver_scheme .NE. 'LxF' ) .AND. ( solver_scheme .NE. 'KT' ) .AND. &
776  ( solver_scheme .NE. 'GFORCE' ) ) THEN
777 
778  WRITE(*,*) 'WARNING: no correct solver scheme selected',solver_scheme
779  WRITE(*,*) 'Chose between: LxF, GFORCE or KT'
780  stop
781 
782  END IF
783 
784  IF ( ( solver_scheme.EQ.'LxF' ) .OR. ( solver_scheme.EQ.'GFORCE' ) ) THEN
785 
786  max_cfl = 1.0
787 
788  ELSEIF ( solver_scheme .EQ. 'KT' ) THEN
789 
790  IF ( ( comp_cells_x .EQ. 1 ) .OR. ( comp_cells_y .EQ. 1 ) ) THEN
791 
792  max_cfl = 0.50d0
793 
794  ELSE
795 
796  max_cfl = 0.25d0
797 
798  END IF
799 
800  END IF
801 
802 
803  IF ( ( cfl .GT. max_cfl ) .OR. ( cfl .LT. 0.d0 ) ) THEN
804 
805  WRITE(*,*) 'WARNING: wrong value of cfl ',cfl
806  WRITE(*,*) 'Choose a value between 0.0 and ',max_cfl
807  READ(*,*)
808 
809  END IF
810 
811  IF ( verbose_level .GE. 1 ) WRITE(*,*) 'Limiters',limiter(1:n_vars)
812 
813  IF ( ( maxval(limiter(1:n_vars)) .GT. 3 ) .OR. &
814  ( minval(limiter(1:n_vars)) .LT. 0 ) ) THEN
815 
816  WRITE(*,*) 'WARNING: wrong limiter ',limiter(1:n_vars)
817  WRITE(*,*) 'Choose among: none, minmod,superbee,van_leer'
818  stop
819 
820  END IF
821 
822  IF ( reconstr_variables .EQ. 'phys' ) THEN
823 
824  IF ( solid_transport_flag ) THEN
825 
826  IF ( temperature_flag ) THEN
827 
828  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
829  WRITE(*,*) 'h+B,u,v,xs,T'
830 
831  ELSE
832 
833  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
834  WRITE(*,*) 'h+B,u,v,xs'
835 
836  END IF
837 
838  ELSE
839 
840  IF ( temperature_flag ) THEN
841 
842  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
843  WRITE(*,*) 'h+B,u,v,T'
844 
845  ELSE
846 
847  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
848  WRITE(*,*) 'h+B,u,v'
849 
850  END IF
851 
852  END IF
853 
854  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
855 
856  IF ( solid_transport_flag ) THEN
857 
858  IF ( temperature_flag ) THEN
859 
860  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
861  WRITE(*,*) 'h+B,hu,hv,xs,T'
862 
863  ELSE
864 
865  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
866  WRITE(*,*) 'h+B,hu,hv,xs'
867 
868  END IF
869 
870  ELSE
871 
872  IF ( solid_transport_flag ) THEN
873 
874  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
875  WRITE(*,*) 'h+B,hu,hv,T'
876 
877  ELSE
878 
879  WRITE(*,*) 'Linear reconstruction and b. c. applied to variables:'
880  WRITE(*,*) 'h+B,hu,hv'
881 
882  END IF
883 
884  END IF
885 
886  END IF
887 
888 
889  IF ( ( reconstr_coeff .GT. 1.0d0 ) .OR. ( reconstr_coeff .LT. 0.d0 ) ) THEN
890 
891  WRITE(*,*) 'WARNING: wrong value of reconstr_coeff ',reconstr_coeff
892  WRITE(*,*) 'Change the value between 0.0 and 1.0 in the input file'
893  READ(*,*)
894 
895  END IF
896 
897  ! ------- READ boundary_conditions NAMELISTS --------------------------------
898 
899  IF ( comp_cells_x .GT. 1 ) THEN
900 
901  ! --------- West boundary conditions -------------------------------------
902 
903  READ(input_unit,west_boundary_conditions,iostat=ios)
904 
905  IF ( ios .NE. 0 ) THEN
906 
907  WRITE(*,*) 'IOSTAT=',ios
908  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
909  WRITE(*,*) 'Please check the input file'
910  stop
911 
912  ELSE
913 
914  rewind(input_unit)
915 
916  END IF
917 
918  IF ( ( hb_bcw%flag .EQ. -1 ) ) THEN
919 
920  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
921  WRITE(*,*) 'B.C. for h+B not set properly'
922  WRITE(*,*) 'Please check the input file'
923  stop
924 
925  END IF
926 
927  IF ( reconstr_variables .EQ. 'phys' ) THEN
928 
929  IF ( u_bcw%flag .EQ. -1 ) THEN
930 
931  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
932  WRITE(*,*) 'B.C. for velocities not set properly'
933  WRITE(*,*) 'Please check the input file'
934  stop
935 
936  END IF
937 
938 
939  IF ( ( comp_cells_y .GT. 1 ) .AND. ( v_bcw%flag .EQ. -1 ) ) THEN
940 
941  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
942  WRITE(*,*) 'B.C. for velocities not set properly'
943  WRITE(*,*) 'Please check the input file'
944  stop
945 
946  END IF
947 
948  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
949 
950  IF ( hu_bcw%flag .EQ. -1 ) THEN
951 
952  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
953  WRITE(*,*) 'B.C. for velocities not set properly'
954  WRITE(*,*) 'Please check the input file'
955  stop
956 
957  END IF
958 
959  IF ( ( comp_cells_y .GT. 1 ) .AND. ( hv_bcw%flag .EQ. -1 ) ) THEN
960 
961  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
962  WRITE(*,*) 'B.C. for velocities not set properly'
963  WRITE(*,*) 'Please check the input file'
964  stop
965 
966  END IF
967 
968  END IF
969 
970  IF ( ( solid_transport_flag ) .AND. ( xs_bcw%flag .EQ. -1 ) ) THEN
971 
972  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
973  WRITE(*,*) 'B.C. for sediment conentration not set properly'
974  WRITE(*,*) 'Please check the input file'
975  stop
976 
977  END IF
978 
979  IF ( ( temperature_flag ) .AND. ( t_bcw%flag .EQ. -1 ) ) THEN
980 
981  WRITE(*,*) 'ERROR: problem with namelist WEST_BOUNDARY_CONDITIONS'
982  WRITE(*,*) 'B.C. for temperature not set properly'
983  WRITE(*,*) 'Please check the input file'
984  stop
985 
986  END IF
987 
988  ! set the approriate boundary conditions
989  bcw(1) = hb_bcw
990 
991  IF ( reconstr_variables .EQ. 'phys' ) THEN
992 
993  bcw(2) = u_bcw
994  bcw(3) = v_bcw
995 
996  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
997 
998  bcw(2) = hu_bcw
999  bcw(3) = hv_bcw
1000 
1001  END IF
1002 
1003  ! ------------- East boundary conditions --------------------------------
1004 
1005  READ(input_unit,east_boundary_conditions,iostat=ios)
1006 
1007  IF ( ios .NE. 0 ) THEN
1008 
1009  WRITE(*,*) 'IOSTAT=',ios
1010  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1011  WRITE(*,*) 'Please check the input file'
1012  stop
1013 
1014  ELSE
1015 
1016  rewind(input_unit)
1017 
1018  END IF
1019 
1020 
1021  IF ( ( hb_bce%flag .EQ. -1 ) ) THEN
1022 
1023  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1024  WRITE(*,*) 'B.C. for h+B not set properly'
1025  WRITE(*,*) 'Please check the input file'
1026  stop
1027 
1028  END IF
1029 
1030  IF ( reconstr_variables .EQ. 'phys' ) THEN
1031 
1032  IF ( u_bce%flag .EQ. -1 ) THEN
1033 
1034  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1035  WRITE(*,*) 'B.C. for velocities not set properly'
1036  WRITE(*,*) 'Please check the input file'
1037  stop
1038 
1039  END IF
1040 
1041  IF ( ( comp_cells_y .GT. 1 ) .AND. ( v_bce%flag .EQ. -1 ) ) THEN
1042 
1043  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1044  WRITE(*,*) 'B.C. for velocities not set properly'
1045  WRITE(*,*) 'Please check the input file'
1046  stop
1047 
1048  END IF
1049 
1050  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
1051 
1052  IF ( hu_bce%flag .EQ. -1 ) THEN
1053 
1054  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1055  WRITE(*,*) 'B.C. for velocities not set properly'
1056  WRITE(*,*) 'Please check the input file'
1057  stop
1058 
1059  END IF
1060 
1061  IF ( ( comp_cells_y .GT. 1 ) .AND. ( hv_bce%flag .EQ. -1 ) ) THEN
1062 
1063  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1064  WRITE(*,*) 'B.C. for velocities not set properly'
1065  WRITE(*,*) 'Please check the input file'
1066  stop
1067 
1068  END IF
1069 
1070  END IF
1071 
1072  IF ( ( solid_transport_flag ) .AND. ( xs_bce%flag .EQ. -1 ) ) THEN
1073 
1074  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1075  WRITE(*,*) 'B.C. for sediment concentration not set properly'
1076  WRITE(*,*) 'Please check the input file'
1077  stop
1078 
1079  END IF
1080 
1081  IF ( ( temperature_flag ) .AND. ( t_bce%flag .EQ. -1 ) ) THEN
1082 
1083  WRITE(*,*) 'ERROR: problem with namelist EAST_BOUNDARY_CONDITIONS'
1084  WRITE(*,*) 'B.C. for temperature not set properly'
1085  WRITE(*,*) 'Please check the input file'
1086  stop
1087 
1088  END IF
1089 
1090  bce(1) = hb_bce
1091 
1092  IF ( reconstr_variables .EQ. 'phys' ) THEN
1093 
1094  bce(2) = u_bce
1095  bce(3) = v_bce
1096 
1097  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
1098 
1099  bce(2) = hu_bce
1100  bce(3) = hv_bce
1101 
1102  END IF
1103 
1104  END IF
1105 
1106  IF ( comp_cells_y .GT. 1 ) THEN
1107 
1108  ! --------------- South boundary conditions ------------------------------
1109 
1110  READ(input_unit,south_boundary_conditions,iostat=ios)
1111 
1112  IF ( ios .NE. 0 ) THEN
1113 
1114  WRITE(*,*) 'IOSTAT=',ios
1115  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1116  WRITE(*,*) 'Please check the input file'
1117  stop
1118 
1119  ELSE
1120 
1121  rewind(input_unit)
1122 
1123  END IF
1124 
1125 
1126  IF ( ( hb_bcs%flag .EQ. -1 ) ) THEN
1127 
1128  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1129  WRITE(*,*) 'B.C. for h+B not set properly'
1130  WRITE(*,*) 'Please check the input file'
1131  stop
1132 
1133  END IF
1134 
1135  IF ( reconstr_variables .EQ. 'phys' ) THEN
1136 
1137  IF ( ( comp_cells_x .GT. 1 ) .AND. ( u_bcs%flag .EQ. -1 ) ) THEN
1138 
1139  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1140  WRITE(*,*) 'B.C. for velocities not set properly'
1141  WRITE(*,*) 'Please check the input file'
1142  stop
1143 
1144  END IF
1145 
1146  IF ( v_bcs%flag .EQ. -1 ) THEN
1147 
1148  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1149  WRITE(*,*) 'B.C. for velocities not set properly'
1150  WRITE(*,*) 'Please check the input file'
1151  stop
1152 
1153  END IF
1154 
1155  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
1156 
1157  IF ( ( comp_cells_x .GT. 1 ) .AND. ( hu_bcs%flag .EQ. -1 ) ) THEN
1158 
1159  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1160  WRITE(*,*) 'B.C. for velocities not set properly'
1161  WRITE(*,*) 'Please check the input file'
1162  stop
1163 
1164  END IF
1165 
1166  IF ( hv_bcs%flag .EQ. -1 ) THEN
1167 
1168  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1169  WRITE(*,*) 'B.C. for velocities not set properly'
1170  WRITE(*,*) 'Please check the input file'
1171  stop
1172 
1173  END IF
1174 
1175  END IF
1176 
1177  IF ( ( solid_transport_flag ) .AND. ( xs_bcs%flag .EQ. -1 ) ) THEN
1178 
1179  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1180  WRITE(*,*) 'B.C. for sediment concentrations not set properly'
1181  WRITE(*,*) 'Please check the input file'
1182  stop
1183 
1184  END IF
1185 
1186  IF ( ( temperature_flag ) .AND. ( t_bcs%flag .EQ. -1 ) ) THEN
1187 
1188  WRITE(*,*) 'ERROR: problem with namelist SOUTH_BOUNDARY_CONDITIONS'
1189  WRITE(*,*) 'B.C. for temperature not set properly'
1190  WRITE(*,*) 'Please check the input file'
1191  stop
1192 
1193  END IF
1194 
1195 
1196  bcs(1) = hb_bcs
1197 
1198  IF ( reconstr_variables .EQ. 'phys' ) THEN
1199 
1200  bcs(2) = u_bcs
1201  bcs(3) = v_bcs
1202 
1203  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
1204 
1205  bcs(2) = hu_bcs
1206  bcs(3) = hv_bcs
1207 
1208  END IF
1209 
1210  ! ---------------- North boundary conditions ----------------------------
1211 
1212  READ(input_unit,north_boundary_conditions,iostat=ios)
1213 
1214  IF ( ios .NE. 0 ) THEN
1215 
1216  WRITE(*,*) 'IOSTAT=',ios
1217  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1218  WRITE(*,*) 'Please check the input file'
1219  stop
1220 
1221  ELSE
1222 
1223  rewind(input_unit)
1224 
1225  END IF
1226 
1227  IF ( ( hb_bcn%flag .EQ. -1 ) ) THEN
1228 
1229  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1230  WRITE(*,*) 'B.C. for h+B not set properly'
1231  WRITE(*,*) 'Please check the input file'
1232  stop
1233 
1234  END IF
1235 
1236  IF ( reconstr_variables .EQ. 'phys' ) THEN
1237 
1238  IF ( ( comp_cells_x .GT. 1 ) .AND. ( u_bcn%flag .EQ. -1 ) ) THEN
1239 
1240  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1241  WRITE(*,*) 'B.C. for velocities not set properly'
1242  WRITE(*,*) 'Please check the input file'
1243  stop
1244 
1245  END IF
1246 
1247  IF ( v_bcn%flag .EQ. -1 ) THEN
1248 
1249  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1250  WRITE(*,*) 'B.C. for velocities not set properly'
1251  WRITE(*,*) 'Please check the input file'
1252  stop
1253 
1254  END IF
1255 
1256  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
1257 
1258  IF ( ( comp_cells_x .GT. 1 ) .AND. ( hu_bcn%flag .EQ. -1 ) ) THEN
1259 
1260  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1261  WRITE(*,*) 'B.C. for velocities not set properly'
1262  WRITE(*,*) 'Please check the input file'
1263  stop
1264 
1265  END IF
1266 
1267  IF ( hv_bcn%flag .EQ. -1 ) THEN
1268 
1269  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1270  WRITE(*,*) 'B.C. for velocities not set properly'
1271  WRITE(*,*) 'Please check the input file'
1272  stop
1273 
1274  END IF
1275 
1276  END IF
1277 
1278  IF (( solid_transport_flag ) .AND. ( xs_bcn%flag .EQ. -1 ) ) THEN
1279 
1280  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1281  WRITE(*,*) 'B.C. for sediment concentrations not set properly'
1282  WRITE(*,*) 'Please check the input file'
1283  stop
1284 
1285  END IF
1286 
1287  IF ( ( temperature_flag ) .AND. ( t_bcn%flag .EQ. -1 ) ) THEN
1288 
1289  WRITE(*,*) 'ERROR: problem with namelist NORTH_BOUNDARY_CONDITIONS'
1290  WRITE(*,*) 'B.C. for temperature not set properly'
1291  WRITE(*,*) 'Please check the input file'
1292  stop
1293 
1294  END IF
1295 
1296 
1297  bcn(1) = hb_bcn
1298 
1299  IF ( reconstr_variables .EQ. 'phys' ) THEN
1300 
1301  bcn(2) = u_bcn
1302  bcn(3) = v_bcn
1303 
1304  ELSEIF ( reconstr_variables .EQ. 'cons' ) THEN
1305 
1306  bcn(2) = hu_bcn
1307  bcn(3) = hv_bcn
1308 
1309  END IF
1310 
1311  END IF
1312 
1313  IF ( solid_transport_flag ) THEN
1314 
1315  bcw(4) = xs_bcw
1316  bce(4) = xs_bce
1317  bcs(4) = xs_bcs
1318  bcn(4) = xs_bcn
1319 
1320  IF ( temperature_flag ) THEN
1321 
1322  bcw(5) = t_bcw
1323  bce(5) = t_bce
1324  bcs(5) = t_bcs
1325  bcn(5) = t_bcn
1326 
1327  END IF
1328 
1329  ELSE
1330 
1331  IF ( temperature_flag ) THEN
1332 
1333  bcw(4) = t_bcw
1334  bce(4) = t_bce
1335  bcs(4) = t_bcs
1336  bcn(4) = t_bcn
1337 
1338  END IF
1339 
1340  END IF
1341 
1342 
1343  ! ------- READ expl_terms_parameters NAMELIST -------------------------------
1344 
1345  READ(input_unit, expl_terms_parameters,iostat=ios)
1346 
1347  IF ( ios .NE. 0 ) THEN
1348 
1349  WRITE(*,*) 'IOSTAT=',ios
1350  WRITE(*,*) 'ERROR: problem with namelist EXPL_TERMS_PARAMETERS'
1351  WRITE(*,*) 'Please check the input file'
1352  stop
1353 
1354  ELSE
1355 
1356  rewind(input_unit)
1357 
1358  END IF
1359 
1360  IF ( grav .EQ. -1.d0 ) THEN
1361 
1362  WRITE(*,*) 'ERROR: problem with namelist EXPL_TERMS_PARAMETERS'
1363  WRITE(*,*) 'R_SOURCE not set properly'
1364  WRITE(*,*) 'Please check the input file'
1365  stop
1366 
1367  END IF
1368 
1369  IF ( source_flag ) THEN
1370 
1371  IF ( r_source .EQ. 0.d0 ) THEN
1372 
1373  WRITE(*,*) 'ERROR: problem with namelist EXPL_TERMS_PARAMETERS'
1374  WRITE(*,*) 'R_SOURCE =',r_source
1375  WRITE(*,*) 'Please check the input file'
1376  stop
1377 
1378  END IF
1379 
1380  IF ( vfr_source .EQ. 0.d0 ) THEN
1381 
1382  WRITE(*,*) 'ERROR: problem with namelist EXPL_TERMS_PARAMETERS'
1383  WRITE(*,*) 'VFR_SOURCE not set properly'
1384  WRITE(*,*) 'Please check the input file'
1385  stop
1386 
1387  END IF
1388 
1389  IF ( temperature_flag .AND. ( r_source .EQ. 0.d0 ) ) THEN
1390 
1391  WRITE(*,*) 'ERROR: problem with namelist EXPL_TERMS_PARAMETERS'
1392  WRITE(*,*) 'TEMPERATURE_FLAG =',temperature_flag,' R_SOURCE =',r_source
1393  WRITE(*,*) 'Please check the input file'
1394  stop
1395 
1396  END IF
1397 
1398  END IF
1399 
1400  ! ------- READ temperature_parameters NAMELIST ------------------------------
1401 
1402  IF ( temperature_flag ) THEN
1403 
1404  READ(input_unit, temperature_parameters,iostat=ios)
1405 
1406  IF ( ios .NE. 0 ) THEN
1407 
1408  WRITE(*,*) 'IOSTAT=',ios
1409  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1410  WRITE(*,*) 'Please check the input file'
1411  stop
1412 
1413  ELSE
1414 
1415  rewind(input_unit)
1416 
1417  END IF
1418 
1419  IF ( exp_area_fract .EQ. -1.d0 ) THEN
1420 
1421  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1422  WRITE(*,*) 'EXP_AREA_FRACT value not properly set'
1423  WRITE(*,*) 'Please check the input file'
1424  stop
1425 
1426  END IF
1427 
1428  IF ( emissivity .EQ. -1.d0 ) THEN
1429 
1430  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1431  WRITE(*,*) 'EMISSIVITY value not properly set'
1432  WRITE(*,*) 'Please check the input file'
1433  stop
1434 
1435  END IF
1436 
1437  IF ( atm_heat_transf_coeff .EQ. -1.d0 ) THEN
1438 
1439  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1440  WRITE(*,*) 'ATM_HEAT_TRANSF_COEFF value not properly set'
1441  WRITE(*,*) 'Please check the input file'
1442  stop
1443 
1444  END IF
1445 
1446  IF ( thermal_conductivity .EQ. -1.d0 ) THEN
1447 
1448  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1449  WRITE(*,*) 'THERMAL CONDUCTIVITY value not properly set'
1450  WRITE(*,*) 'Please check the input file'
1451  stop
1452 
1453  END IF
1454 
1455  IF ( enne .EQ. -1.d0 ) THEN
1456 
1457  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1458  WRITE(*,*) 'ENNE value not properly set'
1459  WRITE(*,*) 'Please check the input file'
1460  stop
1461 
1462  END IF
1463 
1464  IF ( emme .EQ. -1.d0 ) THEN
1465 
1466  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1467  WRITE(*,*) 'EMME value not properly set'
1468  WRITE(*,*) 'Please check the input file'
1469  stop
1470 
1471  END IF
1472 
1473  IF ( t_env .EQ. -1.d0 ) THEN
1474 
1475  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1476  WRITE(*,*) 'T_ENV value not properly set'
1477  WRITE(*,*) 'Please check the input file'
1478  stop
1479 
1480  END IF
1481 
1482  IF ( t_ground .EQ. -1.d0 ) THEN
1483 
1484  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1485  WRITE(*,*) 'T_GROUND value not properly set'
1486  WRITE(*,*) 'Please check the input file'
1487  stop
1488 
1489  END IF
1490 
1491  IF ( c_p .EQ. -1.d0 ) THEN
1492 
1493  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1494  WRITE(*,*) 'C_P value not properly set'
1495  WRITE(*,*) 'Please check the input file'
1496  stop
1497 
1498  END IF
1499 
1500  IF ( rho .EQ. -1.d0 ) THEN
1501 
1502  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1503  WRITE(*,*) 'RHO value not properly set'
1504  WRITE(*,*) 'Please check the input file'
1505  stop
1506 
1507  END IF
1508 
1509 
1510  IF ( emissivity .EQ. 0.d0 ) THEN
1511 
1512  WRITE(*,*) 'No radiative term: emissivity =',emissivity
1513 
1514  END IF
1515 
1516  IF ( atm_heat_transf_coeff .EQ. 0.d0 ) THEN
1517 
1518  WRITE(*,*) 'No convective term: atm_heat_transf_coeff =', &
1519  atm_heat_transf_coeff
1520 
1521  END IF
1522 
1523  IF ( thermal_conductivity .EQ. 0.d0 ) THEN
1524 
1525  WRITE(*,*) 'No conductive term: thermal_conductivity =', &
1526  thermal_conductivity
1527 
1528  END IF
1529 
1530  IF ( rho .LE. 0.d0 ) THEN
1531 
1532  WRITE(*,*) 'ERROR: problem with namelist TEMPERATURE_PARAMETERS'
1533  WRITE(*,*) 'RHO =' , rho
1534  WRITE(*,*) 'Please check the input file'
1535  stop
1536 
1537  END IF
1538 
1539  END IF
1540 
1541 
1542  ! ------- READ rheology_parameters NAMELIST --------------------------------
1543 
1544  IF ( rheology_flag ) THEN
1545 
1546  READ(input_unit, rheology_parameters,iostat=ios)
1547 
1548  IF ( ios .NE. 0 ) THEN
1549 
1550  WRITE(*,*) 'IOSTAT=',ios
1551  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1552  WRITE(*,*) 'Please check the input file'
1553  stop
1554 
1555  ELSE
1556 
1557  rewind(input_unit)
1558 
1559  END IF
1560 
1561  IF ( rheology_model .EQ. 0 ) THEN
1562 
1563  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1564  WRITE(*,*) 'RHEOLOGY_FLAG' , rheology_flag , 'RHEOLOGY_MODEL =' , &
1565  rheology_model
1566  WRITE(*,*) 'Please check the input file'
1567  stop
1568 
1569  ELSEIF ( rheology_model .EQ. 1 ) THEN
1570 
1571  IF ( ( mu .EQ. -1.d0 ) .AND. ( xi .EQ. -1.d0 ) ) THEN
1572 
1573  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1574  WRITE(*,*) 'RHEOLOGY_MODEL =' , rheology_model
1575  WRITE(*,*) 'MU =' , mu ,' XI =' , xi
1576  WRITE(*,*) 'Please check the input file'
1577  stop
1578 
1579  END IF
1580 
1581  IF ( ( t_ref .NE. -1.d0 ) .OR. ( nu_ref .NE. -1.d0 ) .OR. &
1582  ( visc_par .NE. -1.d0 ) .OR. ( tau .NE. -1.d0 ) ) THEN
1583 
1584  WRITE(*,*) 'WARNING: parameters not used in RHEOLOGY_PARAMETERS'
1585  IF ( t_ref .NE. -1.d0 ) WRITE(*,*) 'T_ref =',t_ref
1586  IF ( nu_ref .NE. -1.d0 ) WRITE(*,*) 'nu_ref =',nu_ref
1587  IF ( visc_par .NE. -1.d0 ) WRITE(*,*) 'visc_par =',visc_par
1588  IF ( tau .NE. -1.d0 ) WRITE(*,*) 'tau =',tau
1589  WRITE(*,*) 'Press ENTER to continue'
1590  READ(*,*)
1591 
1592  END IF
1593 
1594  ELSEIF ( rheology_model .EQ. 2 ) THEN
1595 
1596  IF ( tau .EQ. -1.d0 ) THEN
1597 
1598  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1599  WRITE(*,*) 'RHEOLOGY_MODEL =' , rheology_model
1600  WRITE(*,*) 'TAU =' , tau
1601  WRITE(*,*) 'Please check the input file'
1602  stop
1603 
1604  END IF
1605 
1606  IF ( ( t_ref .NE. -1.d0 ) .OR. ( nu_ref .NE. -1.d0 ) .OR. &
1607  ( visc_par .NE. -1.d0 ) .OR. ( mu .NE. -1.d0 ) .OR. &
1608  ( xi .NE. -1.d0 ) ) THEN
1609 
1610  WRITE(*,*) 'WARNING: parameters not used in RHEOLOGY_PARAMETERS'
1611  IF ( t_ref .NE. -1.d0 ) WRITE(*,*) 'T_ref =',t_ref
1612  IF ( nu_ref .NE. -1.d0 ) WRITE(*,*) 'nu_ref =',nu_ref
1613  IF ( visc_par .NE. -1.d0 ) WRITE(*,*) 'visc_par =',visc_par
1614  IF ( mu .NE. -1.d0 ) WRITE(*,*) 'mu =',mu
1615  IF ( xi .NE. -1.d0 ) WRITE(*,*) 'xi =',xi
1616  WRITE(*,*) 'Press ENTER to continue'
1617  READ(*,*)
1618 
1619 
1620  END IF
1621 
1622  ELSEIF ( rheology_model .EQ. 3 ) THEN
1623 
1624  IF ( nu_ref .EQ. -1.d0 ) THEN
1625 
1626  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1627  WRITE(*,*) 'NU_REF =' , nu_ref
1628  WRITE(*,*) 'Please check the input file'
1629  stop
1630 
1631  END IF
1632 
1633  IF ( visc_par .EQ. -1.d0 ) THEN
1634 
1635  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1636  WRITE(*,*) 'VISC_PAR =' , visc_par
1637  WRITE(*,*) 'Please check the input file'
1638  stop
1639 
1640  ELSEIF ( visc_par .EQ. 0.d0 ) THEN
1641 
1642  WRITE(*,*) 'WARNING: temperature and momentum uncoupled'
1643  WRITE(*,*) 'VISC_PAR =' , visc_par
1644  WRITE(*,*) 'Press ENTER to continue'
1645  READ(*,*)
1646 
1647  ELSE
1648 
1649  IF ( t_ref .EQ. -1.d0 ) THEN
1650 
1651  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1652  WRITE(*,*) 'T_REF =' , t_ref
1653  WRITE(*,*) 'Please check the input file'
1654  stop
1655 
1656  END IF
1657 
1658  END IF
1659 
1660  IF ( ( mu .NE. -1.d0 ) .OR. ( xi .NE. -1.d0 ) .OR. ( tau .NE. -1.d0 ) ) &
1661  THEN
1662 
1663  WRITE(*,*) 'WARNING: parameters not used in RHEOLOGY_PARAMETERS'
1664  IF ( mu .NE. -1.d0 ) WRITE(*,*) 'mu =',mu
1665  IF ( xi .NE. -1.d0 ) WRITE(*,*) 'xi =',xi
1666  IF ( tau .NE. -1.d0 ) WRITE(*,*) 'tau =',tau
1667  WRITE(*,*) 'Press ENTER to continue'
1668  READ(*,*)
1669 
1670  END IF
1671 
1672  ELSEIF ( rheology_model .EQ. 4 ) THEN
1673 
1674  IF ( .NOT. solid_transport_flag ) THEN
1675 
1676  IF ( sed_vol_perc .EQ. -1.d0 ) THEN
1677 
1678  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1679  WRITE(*,*) 'RHEOLOGY_MODEL =' , rheology_model
1680  WRITE(*,*) 'SED_VOL_PERC = ' , sed_vol_perc
1681  stop
1682 
1683  END IF
1684 
1685  END IF
1686 
1687  IF ( alpha2 .EQ. -1.d0 ) THEN
1688 
1689  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1690  WRITE(*,*) 'ALPHA2 =' , alpha2
1691  WRITE(*,*) 'Please check the input file'
1692  stop
1693 
1694  END IF
1695 
1696  IF ( beta2 .EQ. -1.d0 ) THEN
1697 
1698  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1699  WRITE(*,*) 'BETA2 =' , beta2
1700  WRITE(*,*) 'Please check the input file'
1701  stop
1702 
1703  END IF
1704 
1705  IF ( alpha1 .EQ. -1.d0 ) THEN
1706 
1707  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1708  WRITE(*,*) 'ALPHA1 =' , alpha1
1709  WRITE(*,*) 'Please check the input file'
1710  stop
1711 
1712  END IF
1713 
1714  IF ( beta1 .EQ. -1.d0 ) THEN
1715 
1716  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1717  WRITE(*,*) 'BETA1 =' , beta1
1718  WRITE(*,*) 'Please check the input file'
1719  stop
1720 
1721  END IF
1722 
1723  IF ( kappa .EQ. -1.d0 ) THEN
1724 
1725  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1726  WRITE(*,*) 'KAPPA =' , kappa
1727  WRITE(*,*) 'Please check the input file'
1728  stop
1729 
1730  END IF
1731 
1732  IF ( n_td .EQ. -1.d0 ) THEN
1733 
1734  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1735  WRITE(*,*) 'N_TD =' , n_td
1736  WRITE(*,*) 'Please check the input file'
1737  stop
1738 
1739  END IF
1740 
1741  IF ( gamma_w .EQ. -1.d0 ) THEN
1742 
1743  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1744  WRITE(*,*) 'GAMMA_W =' , gamma_w
1745  WRITE(*,*) 'Please check the input file'
1746  stop
1747 
1748  END IF
1749 
1750  IF ( gamma_s .EQ. -1.d0 ) THEN
1751 
1752  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1753  WRITE(*,*) 'GAMMA_S =' , gamma_s
1754  WRITE(*,*) 'Please check the input file'
1755  stop
1756 
1757  END IF
1758 
1759  ELSEIF ( rheology_model .EQ. 5 ) THEN
1760 
1761  WRITE(*,*) 'RHEOLOGY_MODEL =' , rheology_model
1762  WRITE(*,*) 'Kurganov & Petrova Example 5'
1763 
1764  ELSE
1765 
1766  WRITE(*,*) 'ERROR: problem with namelist RHEOLOGY_PARAMETERS'
1767  WRITE(*,*) 'RHEOLOGY_MODEL =' , rheology_model
1768  WRITE(*,*) 'Please check the input file'
1769  stop
1770 
1771  END IF
1772 
1773 
1774  END IF
1775 
1776  ! ---------------------------------------------------------------------------
1777 
1778  IF ( verbose_level .GE. 1 ) WRITE(*,*)
1779 
1780  ! read topography from .inp (recommended for simple batimetries)
1781  IF ( .NOT.topography_demfile ) THEN
1782 
1783  tend1 = .false.
1784 
1785  WRITE(*,*) 'Searching for topography_profile'
1786 
1787  topography_profile_search: DO
1788 
1789  READ(input_unit,*, end = 200 ) card
1790 
1791  IF( trim(card) == 'TOPOGRAPHY_PROFILE' ) THEN
1792 
1793  EXIT topography_profile_search
1794 
1795  END IF
1796 
1797  END DO topography_profile_search
1798 
1799 
1800  READ(input_unit,*) n_topography_profile_x
1801 
1802  IF ( verbose_level .GE. 1 ) WRITE(*,*) 'n_topography_profile_x' , &
1803  n_topography_profile_x
1804 
1805  READ(input_unit,*) n_topography_profile_y
1806 
1807  IF ( verbose_level .GE. 1 ) WRITE(*,*) 'n_topography_profile_y' , &
1808  n_topography_profile_y
1809 
1810  ALLOCATE( topography_profile( 3 , n_topography_profile_x , &
1811  n_topography_profile_y ) )
1812 
1813  DO j = 1, n_topography_profile_x
1814 
1815  DO k = 1, n_topography_profile_y
1816 
1817  READ(input_unit,*) topography_profile(1:3,j,k)
1818 
1819  IF ( verbose_level.GE.1 ) WRITE(*,*) j,k,topography_profile(1:3,j,k)
1820 
1821  ENDDO
1822 
1823  END DO
1824 
1825  GOTO 210
1826 200 tend1 = .true.
1827 210 CONTINUE
1828 
1829 
1830  ! read topography from dem file (recommended for complex batimetries)
1831  ELSE
1832 
1833  WRITE(*,*) 'Searching for DEM file'
1834 
1835  INQUIRE(file='topography_dem.asc',exist=lexist)
1836 
1837  IF(lexist)THEN
1838 
1839  OPEN(2001, file='topography_dem.asc', status='old', action='read')
1840 
1841  ELSE
1842 
1843  WRITE(*,*) 'no dem file'
1844  stop
1845 
1846  ENDIF
1847 
1848  READ(2001,*) chara, ncols
1849  READ(2001,*) chara, nrows
1850  READ(2001,*) chara, xllcorner
1851  READ(2001,*) chara, yllcorner
1852  READ(2001,*) chara, cellsize
1853  READ(2001,*) chara, nodata_value
1854 
1855  ! The values read from the DEM files are associated to the center of the
1856  ! pixels. x0 is the left margin of the computational domain and has to be
1857  ! greater than the center of the first pixel.
1858  IF ( x0 .LT. xllcorner + 0.5d0 * cellsize ) THEN
1859 
1860  WRITE(*,*) 'Computational domain problem'
1861  WRITE(*,*) 'x0 < xllcorner+0.5*cellsize',x0,xllcorner+0.5d0*cellsize
1862  stop
1863 
1864  END IF
1865 
1866  ! The right margin of the computational domain should be smaller then the
1867  ! center of the last pixel
1868  IF ( x0 + ( comp_cells_x ) * cell_size .GT. &
1869  xllcorner + ( 0.5d0 + ncols ) * cellsize ) THEN
1870 
1871  WRITE(*,*) 'Computational domain problem'
1872  WRITE(*,*) 'right edge > xllcorner+ncols*cellsize', &
1873  x0+comp_cells_x*cell_size , xllcorner+(0.5d0+ncols)*cellsize
1874  stop
1875 
1876  END IF
1877 
1878  IF ( y0 .LT. yllcorner+0.5d0*cellsize ) THEN
1879 
1880  WRITE(*,*) 'Computational domain problem'
1881  WRITE(*,*) 'y0 < yllcorner+0.5*cellsize',y0,yllcorner+0.5d0*cellsize
1882  stop
1883 
1884  END IF
1885 
1886  IF ( abs( ( y0 + comp_cells_y * cell_size ) - ( yllcorner + 0.5d0 + &
1887  nrows * cellsize ) ) .LT. 1.d-10 ) THEN
1888 
1889  WRITE(*,*) 'Computational domain problem'
1890  WRITE(*,*) 'top edge > yllcorner+nrows*cellsize', &
1891  y0+comp_cells_y*cell_size , yllcorner+(0.5d0+nrows)*cellsize
1892  stop
1893 
1894  END IF
1895 
1896 
1897  WRITE(*,*) 'Reading DEM file'
1898  WRITE(*,*) 'ncols',ncols
1899  WRITE(*,*) 'nrows',nrows
1900 
1901  n_topography_profile_x = ncols
1902 
1903  n_topography_profile_y = nrows
1904 
1905  ALLOCATE( topography_profile( 3 , n_topography_profile_x , &
1906  n_topography_profile_y) )
1907 
1908  DO j=1,n_topography_profile_x
1909 
1910  topography_profile(1,j,:) = xllcorner + ( j - 0.5d0 ) * cellsize
1911 
1912  ENDDO
1913 
1914  DO k=1,n_topography_profile_y
1915 
1916  topography_profile(2,:,k) = yllcorner + ( k - 0.5d0 ) * cellsize
1917 
1918  ENDDO
1919 
1920  ! Read topography values (starting at the upper-left corner)
1921 
1922  DO k=1,n_topography_profile_y
1923 
1924  WRITE(*,fmt="(A1,A,t21,F6.2,A)",advance="NO") achar(13), &
1925  & " Percent Complete: " , &
1926  ( REAL(k) / REAL(n_topography_profile_y))*100.0, "%"
1927 
1928  READ(2001,*) topography_profile(3,:,n_topography_profile_y-k+1)
1929 
1930  ENDDO
1931 
1932  WRITE(*,*) ''
1933 
1934 
1935  CLOSE(2001)
1936 
1937  ENDIF
1938 
1939  ! ------- READ runout_parameters NAMELIST --------------------------------
1940 
1941  IF ( output_runout_flag ) THEN
1942 
1943  READ(input_unit, runout_parameters,iostat=ios)
1944 
1945  IF ( ios .NE. 0 ) THEN
1946 
1947  WRITE(*,*) 'IOSTAT=',ios
1948  WRITE(*,*) 'ERROR: problem with namelist RUNOUT_PARAMETERS'
1949  WRITE(*,*) 'Please check the input file'
1950  stop
1951 
1952  ELSE
1953 
1954  rewind(input_unit)
1955 
1956  END IF
1957 
1958  IF ( ( x0_runout .EQ. -1.d0 ) .AND. ( y0_runout .EQ. -1.d0 ) ) THEN
1959 
1960  WRITE(*,*) 'Runout reference location not defined'
1961 
1962  ELSE
1963 
1964  IF ( x0_runout .LT. x0 ) THEN
1965 
1966  WRITE(*,*) 'Computational domain problem'
1967  WRITE(*,*) 'x0_runout < x0',x0,x0_runout
1968  stop
1969 
1970  END IF
1971 
1972  IF ( x0 .GT. x0+comp_cells_x*cell_size ) THEN
1973 
1974  WRITE(*,*) 'Computational domain problem'
1975  WRITE(*,*) 'x0_runout > x0+comp_cells_x*cell_size' , x0 , &
1976  x0_runout+comp_cells_x*cell_size
1977  stop
1978 
1979  END IF
1980 
1981  IF ( y0_runout .LT. y0 ) THEN
1982 
1983  WRITE(*,*) 'Computational domain problem'
1984  WRITE(*,*) 'y0_runout < y0',y0,y0_runout
1985  stop
1986 
1987  END IF
1988 
1989  IF ( y0 .GT. y0+comp_cells_y*cell_size ) THEN
1990 
1991  WRITE(*,*) 'Computational domain problem'
1992  WRITE(*,*) 'y0_runout > y0+comp_cells_y*cell_size' , y0 , &
1993  y0_runout+comp_cells_y*cell_size
1994  stop
1995 
1996  END IF
1997 
1998  END IF
1999 
2000  runout_file = trim(run_name)//'_runout'//'.txt'
2001 
2002  OPEN(runout_unit,file=runout_file,status='unknown',form='formatted')
2003 
2004  END IF
2005 
2006 
2007  !------ search for check points --------------------------------------------
2008 
2009  tend1 = .false.
2010 
2011  WRITE(*,*) 'Searching for topography_profile'
2012 
2013  n_probes = 0
2014 
2015  probes_search: DO
2016 
2017  READ(input_unit,*, end = 300 ) card
2018 
2019  IF( trim(card) == 'PROBES_COORDS' ) THEN
2020 
2021  EXIT probes_search
2022 
2023  END IF
2024 
2025  END DO probes_search
2026 
2027 
2028  READ(input_unit,*) n_probes
2029 
2030  WRITE(*,*) 'n_probes ',n_probes
2031 
2032  ALLOCATE( probes_coords( 2 , n_probes ) )
2033 
2034  DO k = 1, n_probes
2035 
2036  READ(input_unit,*) probes_coords( 1:2 , k )
2037 
2038  IF ( verbose_level.GE.0 ) WRITE(*,*) k , probes_coords( 1:2 , k )
2039 
2040  END DO
2041 
2042  GOTO 310
2043 300 tend1 = .true.
2044 310 CONTINUE
2045 
2046 
2047 
2048 
2049  ! ----- end search for check points -----------------------------------------
2050 
2051 
2052  CLOSE( input_unit )
2053 
2054  bak_name = trim(run_name)//'.bak'
2055 
2056  OPEN(backup_unit,file=bak_name,status='unknown')
2057 
2058  WRITE(backup_unit, run_parameters )
2059 
2060  IF ( restart ) THEN
2061 
2062  WRITE(backup_unit,restart_parameters)
2063 
2064  ELSE
2065 
2066  WRITE(backup_unit,newrun_parameters)
2067 
2068  IF ( riemann_flag) THEN
2069 
2070  WRITE(backup_unit,left_state)
2071  WRITE(backup_unit,right_state)
2072 
2073  ELSE
2074 
2075  WRITE(backup_unit,initial_conditions)
2076 
2077  END IF
2078 
2079  END IF
2080 
2081  WRITE(backup_unit, numeric_parameters )
2082 
2083  IF ( comp_cells_x .GT. 1 ) THEN
2084 
2085  WRITE(backup_unit,west_boundary_conditions)
2086  WRITE(backup_unit,east_boundary_conditions)
2087 
2088  END IF
2089 
2090  IF ( comp_cells_y .GT. 1 ) THEN
2091 
2092  WRITE(backup_unit,north_boundary_conditions)
2093  WRITE(backup_unit,south_boundary_conditions)
2094 
2095  END IF
2096 
2097  WRITE(backup_unit, expl_terms_parameters )
2098 
2099  IF ( temperature_flag ) WRITE(backup_unit,temperature_parameters)
2100 
2101  IF ( rheology_flag ) WRITE(backup_unit,rheology_parameters)
2102 
2103  IF ( output_runout_flag ) WRITE(backup_unit, runout_parameters)
2104 
2105 
2106  IF ( .NOT.topography_demfile ) THEN
2107 
2108  WRITE(backup_unit,*) '''TOPOGRAPHY_PROFILE'''
2109  WRITE(backup_unit,*) n_topography_profile_x
2110  WRITE(backup_unit,*) n_topography_profile_y
2111 
2112  DO j = 1, n_topography_profile_x
2113 
2114  DO k = 1, n_topography_profile_y
2115 
2116  WRITE(backup_unit,107) topography_profile(1:3,j,k)
2117 
2118 107 FORMAT(3(1x,e14.7))
2119 
2120  ENDDO
2121 
2122  END DO
2123 
2124  END IF
2125 
2126  IF ( n_probes .GT. 0 ) THEN
2127 
2128  WRITE(backup_unit,*) '''PROBES_COORDS'''
2129  WRITE(backup_unit,*) n_probes
2130 
2131  DO k = 1,n_probes
2132 
2133  WRITE(backup_unit,109) probes_coords(1:2,k)
2134 
2135 109 FORMAT(2(1x,e14.7))
2136 
2137  END DO
2138 
2139  END IF
2140 
2141  CLOSE(backup_unit)
2142 
2143  END SUBROUTINE read_param
2144 
2145  !******************************************************************************
2147  !
2152  !
2156  !
2157  !******************************************************************************
2158 
2159  SUBROUTINE update_param
2160 
2161  IMPLICIT none
2162 
2163  INTEGER :: ios
2164 
2165  CHARACTER(LEN=40) :: run_name_org
2166  LOGICAL :: restart_org
2167  LOGICAL :: topography_demfile_org
2168  REAL*8 :: t_start_org
2169  REAL*8 :: t_end_org
2170  REAL*8 :: dt_output_org
2171  LOGICAL :: output_cons_flag_org
2172  LOGICAL :: output_phys_flag_org
2173  LOGICAL :: output_esri_flag_org
2174  INTEGER :: verbose_level_org
2175 
2176  run_name_org = run_name
2177  restart_org = restart
2178  topography_demfile_org = topography_demfile
2179  t_start_org = t_start
2180  t_end_org = t_end
2181  dt_output_org = dt_output
2182  output_cons_flag_org = output_cons_flag
2183  output_phys_flag_org = output_phys_flag
2184  output_esri_flag_org = output_esri_flag
2185  verbose_level_org = verbose_level
2186 
2187 
2188  OPEN(input_unit,file=input_file,status='old')
2189 
2190  ! ------- READ run_parameters NAMELIST -----------------------------------
2191  READ(input_unit, run_parameters,iostat=ios )
2192 
2193  IF ( ios .NE. 0 ) THEN
2194 
2195  WRITE(*,*) 'IOSTAT=',ios
2196  WRITE(*,*) 'ERROR: problem with namelist RUN_PARAMETERS'
2197  WRITE(*,*) 'Please check the input file'
2198  stop
2199 
2200  ELSE
2201 
2202  rewind(input_unit)
2203 
2204  END IF
2205 
2206  CLOSE(input_unit)
2207 
2208  IF ( t_end_org .NE. t_end ) THEN
2209 
2210  WRITE(*,*) 'Modified input file: t_end =',t_end
2211 
2212  END IF
2213 
2214  IF ( dt_output_org .NE. dt_output ) THEN
2215 
2216  WRITE(*,*) 'Modified input file: dt_output =',dt_output
2217 
2218  END IF
2219 
2220  IF ( output_cons_flag_org .NEQV. output_cons_flag ) THEN
2221 
2222  WRITE(*,*) 'Modified input file: output_cons_flag =',output_cons_flag
2223 
2224  END IF
2225 
2226  IF ( output_phys_flag_org .NEQV. output_phys_flag ) THEN
2227 
2228  WRITE(*,*) 'Modified input file: output_phys_flag =',output_phys_flag
2229 
2230  END IF
2231 
2232  IF ( output_esri_flag_org .NEQV. output_esri_flag ) THEN
2233 
2234  WRITE(*,*) 'Modified input file: output_esri_flag =',output_esri_flag
2235 
2236  END IF
2237 
2238  IF ( verbose_level_org .NE. verbose_level ) THEN
2239 
2240  WRITE(*,*) 'Modified input file: verbose_level =',verbose_level
2241 
2242  END IF
2243 
2244  run_name_org = run_name
2245  restart_org = restart
2246  topography_demfile_org = topography_demfile
2247  t_start_org = t_start
2248 
2249 
2250  END SUBROUTINE update_param
2251 
2252 
2253  !******************************************************************************
2255  !
2258  !
2262  !
2263  !******************************************************************************
2264 
2265  SUBROUTINE read_solution
2266 
2267  USE solver_2d, ONLY : allocate_solver_variables
2268 
2269  USE geometry_2d, ONLY : comp_cells_x , x0 , comp_cells_y , y0 , dx , dy
2270  USE geometry_2d, ONLY : b_cent
2271  USE init_2d, ONLY : thickness_init
2272  ! USE parameters_2d, ONLY : n_vars
2273  USE solver_2d, ONLY : q
2274 
2275  IMPLICIT none
2276 
2277  CHARACTER(LEN=15) :: chara
2278 
2279  INTEGER :: j,k
2280 
2281  INTEGER :: dot_idx
2282 
2283  LOGICAL :: lexist
2284 
2285  CHARACTER(LEN=30) :: string
2286 
2287  CHARACTER(LEN=3) :: check_file
2288 
2289  INTEGER :: ncols , nrows , nodata_value
2290 
2291  REAL*8 :: xllcorner , yllcorner , cellsize
2292 
2293  REAL*8 :: xj , yk
2294 
2295  REAL*8 :: sed_vol_fract_init
2296 
2297  INQUIRE (file=restart_file,exist=lexist)
2298 
2299  IF (lexist .EQV. .false.) THEN
2300 
2301  WRITE(*,*) 'Restart: ',trim(restart_file),' not found'
2302  stop
2303 
2304  ELSE
2305 
2306  OPEN(restart_unit,file=restart_file,status='old')
2307 
2308  WRITE(*,*) 'Restart: ',trim(restart_file), ' found'
2309 
2310  END IF
2311 
2312  dot_idx = scan(restart_file, ".", .true.)
2313 
2314  check_file = restart_file(dot_idx+1:dot_idx+3)
2315 
2316  IF ( check_file .EQ. 'asc' ) THEN
2317 
2318  READ(restart_unit,*) chara, ncols
2319  READ(restart_unit,*) chara, nrows
2320  READ(restart_unit,*) chara, xllcorner
2321  READ(restart_unit,*) chara, yllcorner
2322  READ(restart_unit,*) chara, cellsize
2323  READ(restart_unit,*) chara, nodata_value
2324 
2325  ALLOCATE( thickness_init(ncols,nrows) )
2326 
2327  IF ( ncols .NE. comp_cells_x ) THEN
2328 
2329  WRITE(*,*) 'ncols not equal to comp_cells_x',ncols,comp_cells_x
2330  stop
2331 
2332  END IF
2333 
2334  IF ( nrows .NE. comp_cells_y ) THEN
2335 
2336  WRITE(*,*) 'nrows not equal to comp_cells_y',nrows,comp_cells_y
2337  stop
2338 
2339  END IF
2340 
2341  IF ( dabs( xllcorner - x0 ) .GT. 1.d-5*cellsize ) THEN
2342 
2343  WRITE(*,*) 'xllcorner not equal to x0', xllcorner , x0
2344  stop
2345 
2346  END IF
2347 
2348  IF ( dabs( yllcorner - y0 ) .GT. 1.d-5*cellsize ) THEN
2349 
2350  WRITE(*,*) 'yllcorner not equal to y0', yllcorner , y0
2351  stop
2352 
2353  END IF
2354 
2355  IF ( cellsize .NE. cell_size ) THEN
2356 
2357  WRITE(*,*) 'cellsize not equal to cell_size', cellsize , cell_size
2358  stop
2359 
2360  END IF
2361 
2362  DO k=1,nrows
2363 
2364  WRITE(*,fmt="(A1,A,t21,F6.2,A)",advance="NO") achar(13), &
2365  & " Percent Complete: ",( REAL(k) / REAL(comp_cells_y))*100.0, "%"
2366 
2367  READ(restart_unit,*) thickness_init(1:ncols,nrows-k+1)
2368 
2369  ENDDO
2370 
2371  WHERE ( thickness_init .EQ. nodata_value )
2372 
2373  thickness_init = 0.d0
2374 
2375  END WHERE
2376 
2377  WRITE(*,*)
2378  WRITE(*,*) 'Total volume =', cellsize**2 * sum(thickness_init)
2379 
2380  q(1,:,:) = b_cent(:,:) + thickness_init(:,:)
2381 
2382  WRITE(*,*) 'Total volume =',cellsize**2 * sum( thickness_init(:,:) )
2383 
2384  q(2,:,:) = 0.d0
2385  q(3,:,:) = 0.d0
2386 
2387  IF ( solid_transport_flag ) THEN
2388 
2389  q(4,:,:) = thickness_init(:,:) * xs_ambient
2390 
2391  WHERE ( thickness_init .GT. 0.d0 )
2392 
2393  q(4,:,:) = thickness_init(:,:) * xs_init
2394 
2395  END WHERE
2396 
2397  WRITE(*,*) 'MAXVAL(q(4,:,:))',maxval(q(4,:,:))
2398  WRITE(*,*) 'xs_ambient,xs_init',xs_ambient,xs_init
2399 
2400 
2401  sed_vol_fract_init = xs_init * gamma_w / ( xs_init * gamma_w + &
2402  ( 1.d0 - xs_init ) * gamma_s )
2403 
2404  WRITE(*,*) 'Total sediment volume =',cellsize**2*sum( thickness_init* &
2405  sed_vol_fract_init )
2406 
2407  IF ( temperature_flag ) THEN
2408 
2409  q(5,:,:) = thickness_init(:,:) * t_init
2410 
2411  WHERE ( thickness_init .GT. 0.d0 )
2412 
2413  q(5,:,:) = thickness_init(:,:) * t_ambient
2414 
2415  END WHERE
2416 
2417  END IF
2418 
2419  ELSE
2420 
2421  IF ( temperature_flag ) THEN
2422 
2423  q(4,:,:) = thickness_init(:,:) * t_init
2424 
2425  WHERE ( thickness_init .GT. 0.d0 )
2426 
2427  q(4,:,:) = thickness_init(:,:) * t_ambient
2428 
2429  END WHERE
2430 
2431  END IF
2432 
2433  END IF
2434 
2435  output_idx = 0
2436 
2437  IF ( verbose_level .GE. 1 ) THEN
2438 
2439  WRITE(*,*) 'Min q(1,:,:) =',minval(q(1,:,:))
2440  WRITE(*,*) 'Max q(1,:,:) =',maxval(q(1,:,:))
2441  WRITE(*,*) 'SUM(q(1,:,:)) =',sum(q(1,:,:))
2442 
2443  DO k=1,nrows
2444 
2445  WRITE(*,*) k,b_cent(:,k)
2446  READ(*,*)
2447 
2448  END DO
2449 
2450  WRITE(*,*) 'SUM(B_cent(:,:)) =',sum(b_cent(:,:))
2451  READ(*,*)
2452 
2453  END IF
2454 
2455  ELSEIF ( check_file .EQ. 'q_2' ) THEN
2456 
2457  DO k=1,comp_cells_y
2458 
2459  DO j=1,comp_cells_x
2460 
2461  IF ( temperature_flag ) THEN
2462 
2463  READ(restart_unit,1004) xj , yk , q(:,j,k)
2464 
2465  ELSE
2466 
2467  READ(restart_unit,1003) xj , yk , q(:,j,k)
2468 
2469  END IF
2470 
2471  IF ( verbose_level .GE. 2 ) THEN
2472 
2473  WRITE(*,*) 'x,y,q,B',xj,yk,q(:,j,k),b_cent(j,k)
2474 
2475  WRITE(*,*) 'h',q(1,j,k)-b_cent(j,k)
2476 
2477  IF ( solid_transport_flag ) THEN
2478 
2479  WRITE(*,*) 'xs', q(4,j,k) / ( q(1,j,k)-b_cent(j,k) )
2480 
2481  IF ( temperature_flag ) THEN
2482 
2483  WRITE(*,*) 'T', q(5,j,k) / ( q(1,j,k)-b_cent(j,k) )
2484 
2485  END IF
2486 
2487  READ(*,*)
2488 
2489  ELSE
2490 
2491  IF ( temperature_flag ) THEN
2492 
2493  WRITE(*,*) 'T', q(4,j,k) / ( q(1,j,k)-b_cent(j,k) )
2494 
2495  END IF
2496 
2497  READ(*,*)
2498 
2499  END IF
2500 
2501  END IF
2502 
2503  IF ( q(1,j,k) .LE. b_cent(j,k)+1.d-10 ) THEN
2504 
2505  IF ( verbose_level .GE. 2 ) THEN
2506 
2507  WRITE(*,*) 'q(1,j,k),B_cent(j,k)',j,k,q(1,j,k),b_cent(j,k)
2508  READ(*,*)
2509 
2510  END IF
2511 
2512  q(1,j,k) = b_cent(j,k)
2513 
2514  END IF
2515 
2516  ENDDO
2517 
2518  READ(restart_unit,*)
2519 
2520  END DO
2521 
2522  WRITE(*,*) 'Total volume =',dx*dy* sum( q(1,:,:)-b_cent(:,:) )
2523 
2524  IF ( solid_transport_flag ) THEN
2525 
2526  WRITE(*,*) 'Total sediment volume =',dx*dy* sum( q(4,:,:) * &
2527  ( q(1,:,:)-b_cent(:,:) ) )
2528 
2529  END IF
2530 
2531 1003 FORMAT(6e20.12)
2532 1004 FORMAT(7e20.12)
2533 
2534  j = scan(restart_file, '.' , .true. )
2535 
2536  string = trim(restart_file(j-4:j-1))
2537  READ( string,* ) output_idx
2538 
2539  WRITE(*,*)
2540  WRITE(*,*) 'Starting from output index ',output_idx
2541 
2542  ! Set this flag to 0 to not overwrite the initial condition
2543 
2544  ELSE
2545 
2546  WRITE(*,*) 'Restart file not in the right format (*.asc or *)'
2547  stop
2548 
2549  END IF
2550 
2551  CLOSE(restart_unit)
2552 
2553 
2554  END SUBROUTINE read_solution
2555 
2556  !******************************************************************************
2558  !
2563  !
2565  !
2569  !
2570  !******************************************************************************
2571 
2572  SUBROUTINE output_solution(time)
2573 
2574  ! external procedures
2575  USE constitutive_2d, ONLY : qc_to_qp
2576 
2577  ! external variables
2578  USE constitutive_2d, ONLY : h , u , v , xs , t
2579 
2580  USE geometry_2d, ONLY : comp_cells_x , b_cent , comp_cells_y , x_comp, y_comp
2581  ! USE geometry_2d, ONLY : x0 , dx , B_ver , y0 , dy
2582 
2583  USE parameters_2d, ONLY : n_vars
2584  USE parameters_2d, ONLY : t_output , dt_output
2585 
2586  USE solver_2d, ONLY : q
2587 
2588  IMPLICIT none
2589 
2590  REAL*8, INTENT(IN) :: time
2591 
2592  CHARACTER(LEN=4) :: idx_string
2593 
2594  REAL*8 :: qp(n_vars)
2595 
2596  INTEGER :: j,k
2597  INTEGER :: i
2598 
2599  output_idx = output_idx + 1
2600 
2601  idx_string = lettera(output_idx-1)
2602 
2603  IF ( output_cons_flag ) THEN
2604 
2605  output_file_2d = trim(run_name)//'_'//idx_string//'.q_2d'
2606 
2607  WRITE(*,*) 'WRITING ',output_file_2d
2608 
2609  OPEN(output_unit_2d,file=output_file_2d,status='unknown',form='formatted')
2610 
2611  !WRITE(output_unit_2d,1002) x0,dx,comp_cells_x,y0,dy,comp_cells_y,t
2612 
2613  DO k = 1,comp_cells_y
2614 
2615  DO j = 1,comp_cells_x
2616 
2617  DO i = 1,n_vars
2618 
2619  ! Exponents with more than 2 digits cause problems reading
2620  ! into matlab... reset tiny values to zero:
2621  IF ( dabs(q(i,j,k)) .LT. 1d-99) q(i,j,k) = 0.d0
2622 
2623  ENDDO
2624 
2625  IF ( temperature_flag ) THEN
2626 
2627  WRITE(output_unit_2d,1007) x_comp(j), y_comp(k), q(:,j,k)
2628 
2629  ELSE
2630 
2631  WRITE(output_unit_2d,1008) x_comp(j), y_comp(k), q(:,j,k)
2632 
2633  END IF
2634 
2635  ENDDO
2636 
2637  WRITE(output_unit_2d,*) ' '
2638 
2639  END DO
2640 
2641  WRITE(output_unit_2d,*) ' '
2642  WRITE(output_unit_2d,*) ' '
2643 
2644  CLOSE(output_unit_2d)
2645 
2646  END IF
2647 
2648  IF ( output_phys_flag ) THEN
2649 
2650  output_file_2d = trim(run_name)//'_'//idx_string//'.p_2d'
2651 
2652  WRITE(*,*) 'WRITING ',output_file_2d
2653 
2654  OPEN(output_unit_2d,file=output_file_2d,status='unknown',form='formatted')
2655 
2656  DO k = 1,comp_cells_y
2657 
2658  DO j = 1,comp_cells_x
2659 
2660  CALL qc_to_qp(q(:,j,k),b_cent(j,k),qp(:))
2661 
2662  IF ( dabs(REAL(h)) .LT. 1d-99) h = dcmplx(0.d0,0.d0)
2663  IF ( dabs(REAL(u)) .LT. 1d-99) u = dcmplx(0.d0,0.d0)
2664  IF ( dabs(REAL(v)) .LT. 1d-99) v = dcmplx(0.d0,0.d0)
2665 
2666  IF ( solid_transport_flag ) THEN
2667 
2668  IF ( dabs(REAL(xs)) .LT. 1d-99) xs = dcmplx(0.d0,0.d0)
2669 
2670  IF ( temperature_flag ) THEN
2671 
2672  IF ( dabs(REAL(t)) .LT. 1d-99) t = dcmplx(0.d0,0.d0)
2673 
2674  WRITE(output_unit_2d,1010) x_comp(j), y_comp(k), REAL(h), &
2675  REAL(u), REAL(v) , B_cent(j,k) , REAL(h) + B_cent(j,k) , &
2676  REAL(xs) , REAL(t)
2677 
2678  ELSE
2679 
2680  WRITE(output_unit_2d,1009) x_comp(j), y_comp(k), REAL(h), &
2681  REAL(u), REAL(v) , B_cent(j,k) , REAL(h) + B_cent(j,k) , &
2682  REAL(xs)
2683 
2684  END IF
2685 
2686  ELSE
2687 
2688  IF ( temperature_flag ) THEN
2689 
2690  IF ( dabs(REAL(t)) .LT. 1d-99) t = dcmplx(0.d0,0.d0)
2691 
2692  WRITE(output_unit_2d,1012) x_comp(j), y_comp(k), REAL(h), &
2693  REAL(u), REAL(v) , B_cent(j,k) , REAL(h) + B_cent(j,k) , &
2694  REAL(t)
2695 
2696  ELSE
2697 
2698  WRITE(output_unit_2d,1011) x_comp(j), y_comp(k), REAL(h), &
2699  REAL(u), REAL(v) , B_cent(j,k) , REAL(h) + B_cent(j,k)
2700 
2701  END IF
2702 
2703  END IF
2704 
2705  END DO
2706 
2707  WRITE(output_unit_2d,*) ' '
2708 
2709  ENDDO
2710 
2711  WRITE(output_unit_2d,*) ' '
2712  WRITE(output_unit_2d,*) ' '
2713 
2714  CLOSE(output_unit_2d)
2715 
2716  END IF
2717 
2718 1007 FORMAT(7e20.12)
2719 1008 FORMAT(6e20.12)
2720 1009 FORMAT(8e20.12)
2721 1010 FORMAT(9e20.12)
2722 1011 FORMAT(7e20.12)
2723 1012 FORMAT(8e20.12)
2724 
2725  t_output = time + dt_output
2726 
2727  IF ( output_esri_flag ) CALL output_esri(output_idx)
2728 
2729  IF ( n_probes .GT. 0 ) CALL output_probes(output_idx)
2730 
2731  END SUBROUTINE output_solution
2732 
2733  !******************************************************************************
2735  !
2738  !
2740  !
2744  !
2745  !******************************************************************************
2746 
2747  SUBROUTINE output_esri(output_idx)
2748 
2749  USE geometry_2d, ONLY : b_cent , grid_output
2750  ! USE geometry_2d, ONLY : comp_interfaces_x , comp_interfaces_y
2751  USE solver_2d, ONLY : q
2752 
2753  IMPLICIT NONE
2754 
2755  INTEGER, INTENT(IN) :: output_idx
2756 
2757  CHARACTER(LEN=4) :: idx_string
2758  INTEGER :: j
2759 
2760  IF ( output_idx .EQ. 1 ) THEN
2761 
2762  OPEN(dem_esri_unit,file='dem_esri.asc',status='unknown',form='formatted')
2763 
2764  WRITE(dem_esri_unit,'(A,I5)') 'ncols ', comp_cells_x
2765  WRITE(dem_esri_unit,'(A,I5)') 'nrows ', comp_cells_y
2766  WRITE(dem_esri_unit,'(A,F15.3)') 'xllcorner ', x0
2767  WRITE(dem_esri_unit,'(A,F15.3)') 'yllcorner ', y0
2768  WRITE(dem_esri_unit,'(A,F15.3)') 'cellsize ', cell_size
2769  WRITE(dem_esri_unit,'(A,I5)') 'NODATA_value ', -9999
2770 
2771  DO j = comp_cells_y,1,-1
2772 
2773  WRITE(dem_esri_unit,*) b_cent(1:comp_cells_x,j)
2774 
2775  ENDDO
2776 
2777  CLOSE(dem_esri_unit)
2778 
2779  END IF
2780 
2781  idx_string = lettera(output_idx-1)
2782 
2783  !Save thickness
2784  output_esri_file = trim(run_name)//'_'//idx_string//'.asc'
2785 
2786  WRITE(*,*) 'WRITING ',output_esri_file
2787 
2788  OPEN(output_esri_unit,file=output_esri_file,status='unknown',form='formatted')
2789 
2790  grid_output = -9999
2791 
2792  WHERE ( q(1,:,:) - b_cent(:,:) .GE. 1.d-5 )
2793 
2794  grid_output = q(1,:,:) - b_cent(:,:)
2795 
2796  END WHERE
2797 
2798  WRITE(output_esri_unit,'(A,I5)') 'ncols ', comp_cells_x
2799  WRITE(output_esri_unit,'(A,I5)') 'nrows ', comp_cells_y
2800  WRITE(output_esri_unit,'(A,F15.3)') 'xllcorner ', x0
2801  WRITE(output_esri_unit,'(A,F15.3)') 'yllcorner ', y0
2802  WRITE(output_esri_unit,'(A,F15.3)') 'cellsize ', cell_size
2803  WRITE(output_esri_unit,'(A,I5)') 'NODATA_value ', -9999
2804 
2805  DO j = comp_cells_y,1,-1
2806 
2807  WRITE(output_esri_unit,*) grid_output(1:comp_cells_x,j)
2808 
2809  ENDDO
2810 
2811  CLOSE(output_esri_unit)
2812 
2813  IF ( temperature_flag ) THEN
2814 
2815  !Save temperature
2816  output_esri_file = trim(run_name)//'_T_'//idx_string//'.asc'
2817 
2818  WRITE(*,*) 'WRITING ',output_esri_file
2819 
2820  OPEN(output_esri_unit,file=output_esri_file,status='unknown',form='formatted')
2821 
2822  grid_output = -9999
2823 
2824  WHERE ( q(1,:,:) - b_cent(:,:) .GE. 1.d-5 )
2825 
2826  grid_output = q(4,:,:) / ( q(1,:,:) - b_cent(:,:) )
2827 
2828  END WHERE
2829 
2830  WRITE(output_esri_unit,'(A,I5)') 'ncols ', comp_cells_x
2831  WRITE(output_esri_unit,'(A,I5)') 'nrows ', comp_cells_y
2832  WRITE(output_esri_unit,'(A,F15.3)') 'xllcorner ', x0
2833  WRITE(output_esri_unit,'(A,F15.3)') 'yllcorner ', y0
2834  WRITE(output_esri_unit,'(A,F15.3)') 'cellsize ', cell_size
2835  WRITE(output_esri_unit,'(A,I5)') 'NODATA_value ', -9999
2836 
2837  DO j = comp_cells_y,1,-1
2838 
2839  WRITE(output_esri_unit,*) grid_output(1:comp_cells_x,j)
2840 
2841  ENDDO
2842 
2843  CLOSE(output_esri_unit)
2844 
2845  END IF
2846 
2847 
2848  END SUBROUTINE output_esri
2849 
2850  SUBROUTINE close_units
2851 
2852  IMPLICIT NONE
2853 
2854  IF ( output_runout_flag) CLOSE(runout_unit)
2855 
2856  END SUBROUTINE close_units
2857 
2858  !******************************************************************************
2860  !
2863  !
2865  !
2869  !
2870  !******************************************************************************
2871 
2872  CHARACTER*4 FUNCTION lettera(k)
2873  IMPLICIT NONE
2874  CHARACTER ones,tens,hund,thou
2875  !
2876  INTEGER :: k
2877  !
2878  INTEGER :: iten, ione, ihund, ithou
2879  !
2880  ithou=int(k/1000)
2881  ihund=int((k-(ithou*1000))/100)
2882  iten=int((k-(ithou*1000)-(ihund*100))/10)
2883  ione=k-ithou*1000-ihund*100-iten*10
2884  ones=char(ione+48)
2885  tens=char(iten+48)
2886  hund=char(ihund+48)
2887  thou=char(ithou+48)
2888  lettera=thou//hund//tens//ones
2889  !
2890  RETURN
2891  END FUNCTION lettera
2892 
2893  !******************************************************************************
2895  !
2903  !******************************************************************************
2904 
2905  SUBROUTINE output_probes(output_idx)
2906 
2907  USE geometry_2d, ONLY : x_comp , y_comp , b_cent
2908  USE solver_2d, ONLY : q
2909 
2910  USE geometry_2d, ONLY : interp_2d_scalarb
2911 
2912  IMPLICIT NONE
2913 
2914  INTEGER, INTENT(IN) :: output_idx
2915 
2916  CHARACTER(LEN=4) :: idx_string
2917 
2918  REAL*8 :: f2
2919 
2920  INTEGER :: k
2921 
2922  idx_string = lettera(output_idx-1)
2923 
2924  !Save thickness
2925  probes_file = trim(run_name)//'_'//idx_string//'.prb'
2926 
2927  OPEN(probes_unit,file=probes_file,status='unknown',form='formatted')
2928 
2929  DO k=1,n_probes
2930 
2931  CALL interp_2d_scalarb( x_comp , y_comp , q(1,:,:) - b_cent(:,:) , &
2932  probes_coords(1,k) , probes_coords(2,k) , f2 )
2933 
2934  WRITE(probes_unit,'(3e20.12)') probes_coords(1,k) , probes_coords(2,k) , f2
2935 
2936  END DO
2937 
2938  CLOSE(probes_unit)
2939 
2940  END SUBROUTINE output_probes
2941 
2942  !******************************************************************************
2944  !
2953  !******************************************************************************
2954 
2955  SUBROUTINE output_runout(time,stop_flag)
2956 
2957  USE geometry_2d, ONLY : x_comp , y_comp , b_cent
2958  USE parameters_2d, ONLY : t_runout , t_steady
2959  USE solver_2d, ONLY : q
2960 
2961 
2962  IMPLICIT NONE
2963 
2964  REAL*8, INTENT(IN) :: time
2965  LOGICAL, INTENT(INOUT) :: stop_flag
2966 
2967  REAL*8, ALLOCATABLE :: X(:,:), Y(:,:) , dist(:,:)
2968  INTEGER :: sX, sY
2969  INTEGER :: imax(2) , imin(2)
2970  REAL*8 :: max_mom
2971 
2972  sx = size(x_comp)
2973  sy = size(y_comp)
2974 
2975  ALLOCATE( x(sx,sy) , y(sx,sy) , dist(sx,sy) )
2976 
2977  x(:,:) = spread( x_comp, 2, sy )
2978  y(:,:) = spread( y_comp, 1, sx )
2979 
2980  dist(:,:) = 0.d0
2981 
2982  IF ( time .EQ. t_start ) THEN
2983 
2984  ALLOCATE( h_old(sx,sy) )
2985 
2986  h_old(:,:) = q(1,:,:) - b_cent(:,:)
2987 
2988  IF ( ( x0_runout .EQ. -1 ) .AND. ( y0_runout .EQ. -1 ) ) THEN
2989 
2990  WHERE( q(1,:,:) - b_cent(:,:) > 1.d-5 ) dist = b_cent
2991  imin = maxloc( dist )
2992 
2993  x0_runout = x(imin(1),imin(2))
2994  y0_runout = y(imin(1),imin(2))
2995 
2996  WRITE(*,*) 'Runout calculated as linear distance from: (' , &
2997  x0_runout ,',',y0_runout,')'
2998 
2999  WHERE( q(1,:,:) - b_cent(:,:)>1.d-5 ) dist = dsqrt( (x-x0_runout)**2 &
3000  + ( y - y0_runout )**2 )
3001 
3002  imax = maxloc( dist )
3003 
3004  init_runout = dist(imax(1),imax(2))
3005 
3006  END IF
3007 
3008  max_mom = 0.d0
3009 
3010  ELSE
3011 
3012  WHERE( h_old(:,:) > 1.d-5 ) dist = dsqrt( q(2,:,:)**2 + q(3,:,:)**2 )
3013 
3014  max_mom = maxval( dist )
3015 
3016  h_old(:,:) = q(1,:,:) - b_cent(:,:)
3017 
3018  END IF
3019 
3020 
3021  WHERE( q(1,:,:) - b_cent(:,:) > 1.d-5 ) dist = dsqrt( ( x - x0_runout )**2 &
3022  + ( y - y0_runout )**2 )
3023 
3024  imax = maxloc( dist )
3025 
3026  WRITE(runout_unit,'(A,F12.3,A,F12.3,A,F12.3)') 'Time (s) = ',time , &
3027  ' Runout (m) = ',dist(imax(1),imax(2)) - init_runout,' Max mom = ', &
3028  max_mom
3029  CALL flush(runout_unit)
3030 
3031  t_runout = time + dt_runout
3032 
3033  IF ( ( time .GT. t_start ) .AND. ( max_mom .LT. runout_mom_stop ) ) THEN
3034 
3035  WRITE(*,*) 'Steady solution reached, max mom = ', max_mom
3036  stop_flag = .true.
3037 
3038  END IF
3039 
3040  IF ( time .EQ. t_steady ) THEN
3041 
3042  OPEN(dakota_unit,file='dakota.txt',status='unknown',form='formatted')
3043 
3044  WRITE(dakota_unit,*) 'final runout =', dist(imax(1),imax(2)) - init_runout
3045 
3046  CLOSE(dakota_unit)
3047 
3048  END IF
3049 
3050  DEALLOCATE( x , y , dist )
3051 
3052  END SUBROUTINE output_runout
3053 
3054 END MODULE inpout_2d
3055 
constitutive_2d::emissivity
real *8 emissivity
emissivity (eps in Costa & Macedonio, 2005)
Definition: constitutive_2d.f90:72
parameters_2d::max_dt
real *8 max_dt
Largest time step allowed.
Definition: parameters_2d.f90:16
inpout_2d::xllcorner
real *8 xllcorner
Definition: inpout_2d.f90:136
parameters_2d::topography_demfile
logical topography_demfile
Flag for uploading topography from a different file (topography_dem.asc)
Definition: parameters_2d.f90:42
constitutive_2d::qc_to_qp
subroutine qc_to_qp(qc, B, qp)
Conservative to physical variables.
Definition: constitutive_2d.f90:392
inpout_2d::probes_unit
integer, parameter probes_unit
Probes data unit.
Definition: inpout_2d.f90:80
inpout_2d::probes_coords
real *8, dimension(:,:), allocatable probes_coords
Definition: inpout_2d.f90:142
geometry_2d::dy
real *8 dy
Control volumes size.
Definition: geometry_2d.f90:59
parameters_2d::sed_vol_perc
real *8 sed_vol_perc
Definition: parameters_2d.f90:118
parameters_2d::n_rk
integer n_rk
Runge-Kutta order.
Definition: parameters_2d.f90:125
inpout_2d::run_name
character(len=40) run_name
Name of the run.
Definition: inpout_2d.f90:66
parameters_2d::t_steady
real *8 t_steady
end time when reached steady solution
Definition: parameters_2d.f90:153
geometry_2d::x_comp
real *8, dimension(:), allocatable x_comp
Location of the centers (x) of the control volume of the domain.
Definition: geometry_2d.f90:14
constitutive_2d::gamma_s
real *8 gamma_s
Specific weight of sediments.
Definition: constitutive_2d.f90:107
geometry_2d::y0
real *8 y0
Bottom of the physical domain.
Definition: geometry_2d.f90:60
constitutive_2d::nu_ref
real *8 nu_ref
reference kinematic viscosity [m2/s]
Definition: constitutive_2d.f90:51
constitutive_2d::t_ground
real *8 t_ground
temperature of lava-ground interface
Definition: constitutive_2d.f90:78
inpout_2d::x0_runout
real *8 x0_runout
Definition: inpout_2d.f90:148
inpout_2d::output_probes
subroutine output_probes(output_idx)
Write solution at selected points on file.
Definition: inpout_2d.f90:2906
geometry_2d::b_cent
real *8, dimension(:,:), allocatable b_cent
Topography at the centers of the control volumes.
Definition: geometry_2d.f90:35
geometry_2d::comp_cells_x
integer comp_cells_x
Number of control volumes x in the comp. domain.
Definition: geometry_2d.f90:64
constitutive_2d::beta2
real *8 beta2
2nd parameter for yield strenght empirical relationship (O&#39;Brian et al, 1993)
Definition: constitutive_2d.f90:89
init_2d::thickness_init
real *8, dimension(:,:), allocatable thickness_init
Definition: init_2d.f90:17
geometry_2d::y_comp
real *8, dimension(:), allocatable y_comp
Location of the centers (y) of the control volume of the domain.
Definition: geometry_2d.f90:20
parameters_2d::temperature_flag
logical temperature_flag
Flag to choose if we model temperature transport.
Definition: parameters_2d.f90:66
parameters_2d::rheology_flag
logical rheology_flag
Flag to choose if we add the rheology.
Definition: parameters_2d.f90:54
inpout_2d::output_esri_unit
integer, parameter output_esri_unit
Esri Output unit.
Definition: inpout_2d.f90:82
geometry_2d::n_topography_profile_y
integer n_topography_profile_y
Definition: geometry_2d.f90:54
inpout_2d::u_bce
type(bc) u_bce
Definition: inpout_2d.f90:124
inpout_2d::read_solution
subroutine read_solution
Read the solution from the restart unit.
Definition: inpout_2d.f90:2266
inpout_2d::u_bcs
type(bc) u_bcs
Definition: inpout_2d.f90:127
parameters_2d::t_init
real *8 t_init
Initial temperature of the pile of material.
Definition: parameters_2d.f90:107
inpout_2d::xs_bcs
type(bc) xs_bcs
Definition: inpout_2d.f90:127
inpout_2d::v_bcn
type(bc) v_bcn
Definition: inpout_2d.f90:130
inpout_2d::y0_runout
real *8 y0_runout
Definition: inpout_2d.f90:148
init_2d::u_e
real *8 u_e
Right velocity x.
Definition: init_2d.f90:31
constitutive_2d::c_p
real *8 c_p
specific heat [J kg-1 K-1]
Definition: constitutive_2d.f90:60
init_2d::v_w
real *8 v_w
Left velocity y.
Definition: init_2d.f90:26
inpout_2d::v_bcs
type(bc) v_bcs
Definition: inpout_2d.f90:127
inpout_2d::nodata_value
integer nodata_value
Definition: inpout_2d.f90:134
parameters_2d
Parameters.
Definition: parameters_2d.f90:7
geometry_2d::dx
real *8 dx
Control volumes size.
Definition: geometry_2d.f90:56
inpout_2d::hb_bcs
type(bc) hb_bcs
Definition: inpout_2d.f90:127
inpout_2d::output_esri_flag
logical output_esri_flag
Flag to save the output in esri ascii format *.asc.
Definition: inpout_2d.f90:100
parameters_2d::y_release
real *8 y_release
Initial y-coordinate of the pile.
Definition: parameters_2d.f90:94
parameters_2d::t_end
real *8 t_end
end time for the run
Definition: parameters_2d.f90:149
inpout_2d::hb_bce
type(bc) hb_bce
Definition: inpout_2d.f90:124
geometry_2d::comp_cells_y
integer comp_cells_y
Number of control volumes y in the comp. domain.
Definition: geometry_2d.f90:66
solver_2d
Numerical solver.
Definition: solver_2d.f90:12
inpout_2d::ncols
integer ncols
Definition: inpout_2d.f90:134
inpout_2d::restart_file
character(len=40) restart_file
Name of the restart file.
Definition: inpout_2d.f90:70
inpout_2d::hv_bcw
type(bc) hv_bcw
Definition: inpout_2d.f90:121
constitutive_2d::xi
real *8 xi
Definition: constitutive_2d.f90:30
parameters_2d::rheology_model
integer rheology_model
choice of the rheology model
Definition: parameters_2d.f90:60
inpout_2d::probes_file
character(len=40) probes_file
Name of the probes file.
Definition: inpout_2d.f90:71
parameters_2d::bcw
type(bc), dimension(:), allocatable bcw
bcW&flag defines the west boundary condition:
Definition: parameters_2d.f90:172
inpout_2d
Input/Output module.
Definition: inpout_2d.f90:13
inpout_2d::hv_bce
type(bc) hv_bce
Definition: inpout_2d.f90:124
init_2d::riemann_interface
real *8 riemann_interface
Riemann problem interface relative position. It is a value between 0 and 1.
Definition: init_2d.f90:21
init_2d
Initial solution.
Definition: init_2d.f90:8
inpout_2d::update_param
subroutine update_param
Read the input file.
Definition: inpout_2d.f90:2160
inpout_2d::u_bcn
type(bc) u_bcn
Definition: inpout_2d.f90:130
inpout_2d::hv_bcn
type(bc) hv_bcn
Definition: inpout_2d.f90:130
init_2d::u_w
real *8 u_w
Left velocity x.
Definition: init_2d.f90:25
constitutive_2d::thermal_conductivity
real *8 thermal_conductivity
thermal conductivity [W m-1 K-1] (k in Costa & Macedonio, 2005)
Definition: constitutive_2d.f90:81
parameters_2d::velocity_ang_release
real *8 velocity_ang_release
Initial velocity direction (angle in degree): .
Definition: parameters_2d.f90:104
parameters_2d::velocity_mod_release
real *8 velocity_mod_release
Initial velocity module of the pile.
Definition: parameters_2d.f90:97
geometry_2d::n_topography_profile_x
integer n_topography_profile_x
Definition: geometry_2d.f90:54
constitutive_2d::beta1
real *8 beta1
2nd parameter for fluid viscosity empirical relationship (O&#39;Brian et al, 1993)
Definition: constitutive_2d.f90:95
parameters_2d::released_volume
real *8 released_volume
Initial volume of the flow.
Definition: parameters_2d.f90:88
inpout_2d::output_esri_file
character(len=40) output_esri_file
Name of the esri output files.
Definition: inpout_2d.f90:73
parameters_2d::n_vars
integer n_vars
Number of conservative variables.
Definition: parameters_2d.f90:120
constitutive_2d::v
complex *16 v
velocity (y direction)
Definition: constitutive_2d.f90:21
inpout_2d::dakota_unit
integer, parameter dakota_unit
Definition: inpout_2d.f90:85
parameters_2d::cfl
real *8 cfl
Courant-Friedrichs-Lewy parameter.
Definition: parameters_2d.f90:18
parameters_2d::t_ambient
real *8 t_ambient
Ambient temperature.
Definition: parameters_2d.f90:110
constitutive_2d::enne
real *8 enne
thermal boundary layer fraction of total thickness
Definition: constitutive_2d.f90:75
parameters_2d::source_flag
logical source_flag
Flag to choose if there is a source of mass within the domain.
Definition: parameters_2d.f90:78
parameters_2d::vfr_source
real *8 vfr_source
Definition: parameters_2d.f90:83
constitutive_2d::rho
real *8 rho
fluid density [kg/m3]
Definition: constitutive_2d.f90:45
parameters_2d::dt0
real *8 dt0
Initial time step.
Definition: parameters_2d.f90:14
constitutive_2d
Constitutive equations.
Definition: constitutive_2d.f90:4
constitutive_2d::grav
real *8 grav
gravitational acceleration
Definition: constitutive_2d.f90:26
inpout_2d::init_runout
real *8 init_runout
Definition: inpout_2d.f90:148
inpout_2d::t_bcw
type(bc) t_bcw
Definition: inpout_2d.f90:121
constitutive_2d::mu
real *8 mu
drag coefficients (Voellmy-Salm model)
Definition: constitutive_2d.f90:29
constitutive_2d::visc_par
real *8 visc_par
viscosity parameter [K-1] (b in Table 1 Costa & Macedonio, 2005)
Definition: constitutive_2d.f90:54
parameters_2d::x_release
real *8 x_release
Initial x-coordiante of the pile.
Definition: parameters_2d.f90:91
inpout_2d::hb_bcw
type(bc) hb_bcw
Definition: inpout_2d.f90:121
inpout_2d::xs_bcw
type(bc) xs_bcw
Definition: inpout_2d.f90:121
inpout_2d::nrows
integer nrows
Definition: inpout_2d.f90:134
constitutive_2d::n_td
real *8 n_td
Mannings roughness coefficient ( units: T L^(-1/3) )
Definition: constitutive_2d.f90:101
inpout_2d::restart
logical restart
Flag to start a run from a previous output: .
Definition: inpout_2d.f90:94
geometry_2d::x0
real *8 x0
Left of the physical domain.
Definition: geometry_2d.f90:57
inpout_2d::output_runout
subroutine output_runout(time, stop_flag)
Write runout on file.
Definition: inpout_2d.f90:2956
inpout_2d::output_solution
subroutine output_solution(time)
Write the solution on the output unit.
Definition: inpout_2d.f90:2573
inpout_2d::backup_unit
integer, parameter backup_unit
Backup input data unit.
Definition: inpout_2d.f90:77
inpout_2d::runout_mom_stop
real *8 runout_mom_stop
Definition: inpout_2d.f90:148
inpout_2d::lettera
character *4 function lettera(k)
Numeric to String conversion.
Definition: inpout_2d.f90:2873
parameters_2d::solver_scheme
character(len=20) solver_scheme
Finite volume method: .
Definition: parameters_2d.f90:145
inpout_2d::output_phys_flag
logical output_phys_flag
Flag to save the physical variables on file *.p_2d.
Definition: inpout_2d.f90:106
geometry_2d::grid_output
real *8, dimension(:,:), allocatable grid_output
Solution in ascii grid format (ESRI)
Definition: geometry_2d.f90:44
inpout_2d::output_unit_2d
integer, parameter output_unit_2d
Output data 2D unit.
Definition: inpout_2d.f90:81
constitutive_2d::t_env
real *8 t_env
evironment temperature [K]
Definition: constitutive_2d.f90:36
constitutive_2d::alpha2
real *8 alpha2
1st parameter for yield strenght empirical relationship (O&#39;Brian et al, 1993)
Definition: constitutive_2d.f90:86
parameters_2d::x_source
real *8 x_source
Definition: parameters_2d.f90:80
inpout_2d::xs_bcn
type(bc) xs_bcn
Definition: inpout_2d.f90:130
inpout_2d::output_unit
integer, parameter output_unit
Output data unit.
Definition: inpout_2d.f90:78
constitutive_2d::gamma_w
real *8 gamma_w
Specific weight of water.
Definition: constitutive_2d.f90:104
geometry_2d
Grid module.
Definition: geometry_2d.f90:7
inpout_2d::hb_bcn
type(bc) hb_bcn
Definition: inpout_2d.f90:130
inpout_2d::hu_bcs
type(bc) hu_bcs
Definition: inpout_2d.f90:127
geometry_2d::cell_size
real *8 cell_size
Definition: geometry_2d.f90:68
inpout_2d::runout_unit
integer, parameter runout_unit
Definition: inpout_2d.f90:84
inpout_2d::output_cons_flag
logical output_cons_flag
Flag to save the conservative variables on file *.q_2d.
Definition: inpout_2d.f90:112
parameters_2d::theta
real *8 theta
Van Leer limiter parameter.
Definition: parameters_2d.f90:147
parameters_2d::riemann_flag
logical riemann_flag
Flag to choose the sort of problem to solve.
Definition: parameters_2d.f90:48
inpout_2d::output_file
character(len=40) output_file
Name of the output files.
Definition: inpout_2d.f90:69
init_2d::t_w
real *8 t_w
Left temperature.
Definition: init_2d.f90:28
constitutive_2d::atm_heat_transf_coeff
real *8 atm_heat_transf_coeff
atmospheric heat trasnfer coefficient [W m-2 K-1] (lambda in C&M, 2005)
Definition: constitutive_2d.f90:63
init_2d::v_e
real *8 v_e
Right velocity y.
Definition: init_2d.f90:32
parameters_2d::t_output
real *8 t_output
time of the next output
Definition: parameters_2d.f90:150
inpout_2d::u_bcw
type(bc) u_bcw
Definition: inpout_2d.f90:121
inpout_2d::read_param
subroutine read_param
Read the input file.
Definition: inpout_2d.f90:529
inpout_2d::t_bcn
type(bc) t_bcn
Definition: inpout_2d.f90:130
inpout_2d::v_bce
type(bc) v_bce
Definition: inpout_2d.f90:124
constitutive_2d::alpha1
real *8 alpha1
1st parameter for fluid viscosity empirical relationship (O&#39;Brian et al, 1993)
Definition: constitutive_2d.f90:92
parameters_2d::xs_init
real *8 xs_init
Initial sediment concentration in the pile of material.
Definition: parameters_2d.f90:113
geometry_2d::topography_profile
real *8, dimension(:,:,:), allocatable topography_profile
Definition: geometry_2d.f90:52
constitutive_2d::t_ref
real *8 t_ref
reference temperature [K]
Definition: constitutive_2d.f90:48
inpout_2d::input_file
character(len=40) input_file
File with the run parameters.
Definition: inpout_2d.f90:68
parameters_2d::xs_ambient
real *8 xs_ambient
Ambient sediment concentration.
Definition: parameters_2d.f90:116
inpout_2d::hu_bce
type(bc) hu_bce
Definition: inpout_2d.f90:124
parameters_2d::y_source
real *8 y_source
Definition: parameters_2d.f90:81
parameters_2d::reconstr_variables
character(len=20) reconstr_variables
Definition: parameters_2d.f90:22
geometry_2d::interp_2d_scalarb
subroutine interp_2d_scalarb(x1, y1, f1, x2, y2, f2)
Scalar interpolation (2D)
Definition: geometry_2d.f90:507
inpout_2d::xs_bce
type(bc) xs_bce
Definition: inpout_2d.f90:124
inpout_2d::output_esri
subroutine output_esri(output_idx)
Write the thickness in ESRI format.
Definition: inpout_2d.f90:2748
parameters_2d::topography_function_flag
logical topography_function_flag
Flag to choose in which way we upload the topography.
Definition: parameters_2d.f90:36
inpout_2d::dt_runout
real *8 dt_runout
Definition: inpout_2d.f90:144
inpout_2d::hu_bcw
type(bc) hu_bcw
Definition: inpout_2d.f90:121
constitutive_2d::exp_area_fract
real *8 exp_area_fract
fractional area of the exposed inner core (f in C&M, 2005)
Definition: constitutive_2d.f90:66
parameters_2d::limiter
integer, dimension(10) limiter
Limiter for the slope in the linear reconstruction: .
Definition: parameters_2d.f90:138
inpout_2d::hv_bcs
type(bc) hv_bcs
Definition: inpout_2d.f90:127
inpout_2d::output_idx
integer output_idx
Counter for the output files.
Definition: inpout_2d.f90:88
inpout_2d::close_units
subroutine close_units
Definition: inpout_2d.f90:2851
inpout_2d::cellsize
real *8 cellsize
Definition: inpout_2d.f90:136
parameters_2d::bcn
type(bc), dimension(:), allocatable bcn
bcN&flag defines the north boundary condition:
Definition: parameters_2d.f90:202
inpout_2d::t_bce
type(bc) t_bce
Definition: inpout_2d.f90:124
parameters_2d::t_runout
real *8 t_runout
time of the next runout output
Definition: parameters_2d.f90:152
init_2d::t_e
real *8 t_e
Right temperature.
Definition: init_2d.f90:34
inpout_2d::hu_bcn
type(bc) hu_bcn
Definition: inpout_2d.f90:130
parameters_2d::n_eqns
integer n_eqns
Number of equations.
Definition: parameters_2d.f90:121
constitutive_2d::xs
complex *16 xs
sediment concentration
Definition: constitutive_2d.f90:23
inpout_2d::output_file_2d
character(len=40) output_file_2d
Name of the output files.
Definition: inpout_2d.f90:72
inpout_2d::dem_esri_unit
integer, parameter dem_esri_unit
Computational grid Esri fmt unit.
Definition: inpout_2d.f90:83
inpout_2d::bak_name
character(len=40) bak_name
Backup file for the parameters.
Definition: inpout_2d.f90:67
inpout_2d::v_bcw
type(bc) v_bcw
Definition: inpout_2d.f90:121
inpout_2d::init_param
subroutine init_param
Initialization of the variables read from the input file.
Definition: inpout_2d.f90:213
constitutive_2d::h
complex *16 h
height [m]
Definition: constitutive_2d.f90:19
inpout_2d::h_old
real *8, dimension(:,:), allocatable h_old
Definition: inpout_2d.f90:146
init_2d::xs_e
real *8 xs_e
Right sediment concentration.
Definition: init_2d.f90:33
init_2d::xs_w
real *8 xs_w
Left sediment concentration.
Definition: init_2d.f90:27
parameters_2d::reconstr_coeff
real *8 reconstr_coeff
Slope coefficient in the linear reconstruction.
Definition: parameters_2d.f90:24
parameters_2d::bcs
type(bc), dimension(:), allocatable bcs
bcS&flag defines the south boundary condition:
Definition: parameters_2d.f90:192
constitutive_2d::emme
real *8 emme
velocity boundary layer fraction of total thickness
Definition: constitutive_2d.f90:57
inpout_2d::write_first_q
logical write_first_q
Definition: inpout_2d.f90:138
inpout_2d::n_probes
integer n_probes
Definition: inpout_2d.f90:140
inpout_2d::t_bcs
type(bc) t_bcs
Definition: inpout_2d.f90:127
constitutive_2d::tau
real *8 tau
drag coefficients (plastic model)
Definition: constitutive_2d.f90:33
parameters_2d::bce
type(bc), dimension(:), allocatable bce
bcE&flag defines the east boundary condition:
Definition: parameters_2d.f90:182
constitutive_2d::t
complex *16 t
temperature
Definition: constitutive_2d.f90:22
constitutive_2d::u
complex *16 u
velocity (x direction)
Definition: constitutive_2d.f90:20
inpout_2d::yllcorner
real *8 yllcorner
Definition: inpout_2d.f90:136
parameters_2d::dt_output
real *8 dt_output
time interval for the output of the solution
Definition: parameters_2d.f90:151
inpout_2d::input_unit
integer, parameter input_unit
Input data unit.
Definition: inpout_2d.f90:76
constitutive_2d::kappa
real *8 kappa
Empirical resistance parameter.
Definition: constitutive_2d.f90:98
inpout_2d::restart_unit
integer, parameter restart_unit
Restart data unit.
Definition: inpout_2d.f90:79
parameters_2d::solid_transport_flag
logical solid_transport_flag
Flag to choose if we model solid phase transport.
Definition: parameters_2d.f90:72
inpout_2d::runout_file
character(len=40) runout_file
Name of the runout file.
Definition: inpout_2d.f90:74
solver_2d::allocate_solver_variables
subroutine allocate_solver_variables
Memory allocation.
Definition: solver_2d.f90:141
init_2d::hb_w
real *8 hb_w
Left height.
Definition: init_2d.f90:24
init_2d::hb_e
real *8 hb_e
Right height.
Definition: init_2d.f90:30
parameters_2d::interfaces_relaxation
logical interfaces_relaxation
Flag to add the relaxation terms after the linear reconstruction: .
Definition: parameters_2d.f90:30
inpout_2d::output_runout_flag
logical output_runout_flag
Flag to save the max runout at ouput times.
Definition: inpout_2d.f90:118
parameters_2d::t_start
real *8 t_start
initial time for the run
Definition: parameters_2d.f90:148
parameters_2d::r_source
real *8 r_source
Definition: parameters_2d.f90:82
parameters_2d::bc
Definition: parameters_2d.f90:157
parameters_2d::t_source
real *8 t_source
Definition: parameters_2d.f90:85
solver_2d::q
real *8, dimension(:,:,:), allocatable q
Conservative variables.
Definition: solver_2d.f90:33
Generated on Tue Feb 13 2018 12:40:55 for IMEX_SfloW2D by   doxygen 1.8.11