inpout_2d.f90

Go to the documentation of this file.

!********************************************************************************
!
!
!
!********************************************************************************

MODULE inpout_2d

  USE inpout, ONLY : run_name
  
  USE variables, ONLY : verbose_level , gi
  USE variables, ONLY : wp

  USE parameters_2d, ONLY : n_vars

  ! -- Variables for the namelist RUN_PARAMETERS
  USE parameters_2d, ONLY : t_start , t_end , dt_output 

  ! -- Variables for the namelist NEWRUN_PARAMETERS
  USE geometry_2d, ONLY : x0 , y0 , comp_cells_x , comp_cells_y , cell_size , dx , dy

  ! -- Variables for the namelists LEFT/RIGHT_BOUNDARY_CONDITIONS
  USE parameters_2d, ONLY : bc

  ! -- Variables for the namelist NUMERIC_PARAMETERS
  USE parameters_2d, ONLY : rsource_cells , solver_scheme, dt0 , max_dt , cfl, limiter , theta, &
       reconstr_coeff , interfaces_relaxation , n_rk   

  ! -- Variables for the namelist SOURCE_PARAMETERS
  USE parameters_2d, ONLY : x_source , y_source , r_source , vel_source , h_source , time_param
  USE parameters_2d, ONLY : vol_flux_source , u_source , v_source , dr_dz

  ! -- Variables for the namelist ATM_PARAMETERS
  USE constitutive_2d, ONLY : grav , rho_nbl , drho_dz , n

  IMPLICIT NONE

  CHARACTER(LEN=40) :: output_file
  CHARACTER(LEN=40) :: output_file_2d

  INTEGER, PARAMETER :: output_unit = 9
  INTEGER, PARAMETER :: output_unit_2d = 12

  INTEGER :: output_idx 

  LOGICAL :: output_phys_flag

  LOGICAL :: output_cons_flag

  ! -- Variables for the namelists WEST_BOUNDARY_CONDITIONS
  TYPE(bc) :: h_bcw , hu_bcw , hv_bcw

  ! -- Variables for the namelists EAST_BOUNDARY_CONDITIONS
  TYPE(bc) :: h_bce , hu_bce , hv_bce

  ! -- Variables for the namelists SOUTH_BOUNDARY_CONDITIONS
  TYPE(bc) :: h_bcs , hu_bcs , hv_bcs

  ! -- Variables for the namelists NORTH_BOUNDARY_CONDITIONS
  TYPE(bc) :: h_bcn , hu_bcn , hv_bcn

  ! parameters to read a dem file
  INTEGER :: ncols, nrows, nodata_value

  REAL(wp) :: xllcorner, yllcorner, cellsize

  LOGICAL :: write_first_q

  INTEGER :: n_probes

  REAL(wp), ALLOCATABLE :: probes_coords(:,:)

  REAL(wp), ALLOCATABLE :: h_old(:,:)

CONTAINS

  !******************************************************************************
  !
  !
  !
  !******************************************************************************

  SUBROUTINE init_param

    USE parameters_2d, ONLY : n_vars , n_eqns 
    USE parameters_2d, ONLY : limiter
    USE parameters_2d, ONLY : bcw , bce , bcs , bcn

    IMPLICIT none

    LOGICAL :: lexist

    n_vars = 3

    !-- Inizialization of the Variables for the namelist RUN_PARAMETERS
    t_start = 0.0
    output_phys_flag = .true.

    n_vars = 3
    n_eqns = 3

    ALLOCATE( bcw(n_vars) , bce(n_vars) , bcs(n_vars) , bcn(n_vars) )
    
    h_bcw%FLAG=          1
    h_bcw%VALUE=  0.0000000000000000     
    hu_bcw%FLAG=          1
    hu_bcw%VALUE=  0.0000000000000000     
    hv_bcw%FLAG=          1
    hv_bcw%VALUE=  0.0000000000000000 

    bcw(1) = h_bcw
    bcw(2) = hu_bcw 
    bcw(3) = hv_bcw 

    h_bce%FLAG=          1
    h_bce%VALUE=  0.0000000000000000     
    hu_bce%FLAG=          1
    hu_bce%VALUE=  0.0000000000000000     
    hv_bce%FLAG=          1
    hv_bce%VALUE=  0.0000000000000000     

    bce(1) = h_bce 
    bce(2) = hu_bce 
    bce(3) = hv_bce 

    h_bcs%FLAG=          1
    h_bcs%VALUE=  0.0000000000000000     
    hu_bcs%FLAG=          1
    hu_bcs%VALUE=  0.0000000000000000     
    hv_bcs%FLAG=          1
    hv_bcs%VALUE=  0.0000000000000000     

    bcs(1) = h_bcs 
    bcs(2) = hu_bcs 
    bcs(3) = hv_bcs 

    h_bcn%FLAG=          1
    h_bcn%VALUE=  0.0000000000000000     
    hu_bcn%FLAG=          1
    hu_bcn%VALUE=  0.0000000000000000     
    hv_bcn%FLAG=          1
    hv_bcn%VALUE=  0.0000000000000000     

    bcn(1) = h_bcn 
    bcn(2) = hu_bcn 
    bcn(3) = hv_bcn 

    time_param(1) = t_end
    time_param(2) = t_end
    time_param(3) = 0.0_wp
    time_param(4) = t_end
        
    cell_size = r_source / rsource_cells
    comp_cells_x = floor( 40.0_wp * r_source / cell_size )
    comp_cells_y = comp_cells_x
    x0 = x_source - 0.25_wp * comp_cells_x * cell_size
    y0 = y_source - 0.25_wp * comp_cells_y * cell_size
    dx = cell_size
    dy = cell_size
        
    grav = gi
    n = sqrt( - grav / rho_nbl * drho_dz )

    WRITE(*,*) 'x_source =',x_source
    WRITE(*,*) 'y_source =',y_source
    WRITE(*,*) 'r_source =',r_source
    WRITE(*,*) 'grav =',grav
    WRITE(*,*) 'rho_nbl =',rho_nbl
    WRITE(*,*) 'drho_dz =',drho_dz
    
    WRITE(*,*) 'N',n
    
  END SUBROUTINE init_param


  !******************************************************************************
  !
  !
  !
  !
  !******************************************************************************

  SUBROUTINE output_solution(time,steady_state)

    ! external procedures
    USE constitutive_2d, ONLY : qc_to_qp

    USE geometry_2d, ONLY : comp_cells_x , comp_cells_y , x_comp,      &
         y_comp

    USE parameters_2d, ONLY : n_vars
    USE parameters_2d, ONLY : t_output , dt_output 

    USE solver_2d, ONLY : q

    IMPLICIT none

    REAL(wp), INTENT(IN) :: time
    LOGICAL, INTENT(IN) :: steady_state
    
    CHARACTER(LEN=4) :: idx_string

    REAL(wp) :: qp(n_vars+2)

    REAL(wp) :: r_u , r_v , r_h 

    INTEGER :: j,k
    INTEGER :: i
    INTEGER :: i_vars

    IF ( .NOT. steady_state ) THEN
    
       output_idx = output_idx + 1
       
       idx_string = lettera(output_idx-1)
       
       output_file_2d = trim(run_name)//'_'//idx_string//'.p_2d'

    ELSE

       output_file_2d = trim(run_name)//'_'//'std.p_2d'

    END IF

       
    IF ( verbose_level .GE. 0 ) WRITE(*,*) 'WRITING ',output_file_2d

    OPEN(output_unit_2d,file=output_file_2d,status='unknown',form='formatted')

    DO k = 1,comp_cells_y

       DO j = 1,comp_cells_x

          CALL qc_to_qp(q(1:n_vars,j,k) , qp(1:n_vars+2))

          r_h = qp(1)
          r_u = qp(n_vars+1)
          r_v = qp(n_vars+2)

          IF ( abs( r_h ) .LT. 1e-20_wp) r_h = 0.0_wp
          IF ( abs( r_u ) .LT. 1e-20_wp) r_u = 0.0_wp
          IF ( abs( r_v ) .LT. 1e-20_wp) r_v = 0.0_wp

          WRITE(output_unit_2d,1010) x_comp(j), y_comp(k), r_h , r_u , r_v

       END DO

       WRITE(output_unit_2d,*) ' ' 

    ENDDO

    WRITE(output_unit_2d,*) ' '
    WRITE(output_unit_2d,*) ' '

    CLOSE(output_unit_2d)


1010 FORMAT(100es15.7e2)

    t_output = time + dt_output

  END SUBROUTINE output_solution

  

  SUBROUTINE close_units

    IMPLICIT NONE

  END SUBROUTINE close_units

  !******************************************************************************
  !
  !
  !
  !
  !******************************************************************************

  CHARACTER*4 FUNCTION lettera(k)
    IMPLICIT NONE
    CHARACTER ones,tens,hund,thou
    !
    INTEGER :: k
    !
    INTEGER :: iten, ione, ihund, ithou
    !
    ithou=int(k/1000)
    ihund=int((k-(ithou*1000))/100)
    iten=int((k-(ithou*1000)-(ihund*100))/10)
    ione=k-ithou*1000-ihund*100-iten*10
    ones=char(ione+48)
    tens=char(iten+48)
    hund=char(ihund+48)
    thou=char(ithou+48)
    lettera=thou//hund//tens//ones
    !
    RETURN
  END FUNCTION lettera

END MODULE inpout_2d