module model_zhdanov_m !! Implementation of ZHDANOV model use precision_grillix_m, only : GP, MPI_GP use comm_handler_m, only : comm_handler_t use parallelisation_setup_m, only : is_rank_info_writer use static_data_m, only : equi, & equi_storage => equi_on_cano, & mesh => mesh_cano, & multigrid => multigrid_cano, & map, & polars => polars_cano, & penalisation => parbnd_immersed_cano use variable_m, only : variable_t use multistep_m, only : karniadakis_t use checkpoint_monitor_m, only : checkpoint_monitor_t use snapshot_m, only : snapshot_t, char_arr_t use init_zhdanov_m, only : build_zhdanov_initial_state, & initialise_timesteppers use timestep_zhdanov_m, only : timestep_zhdanov use params_zhdanov_general_m, only : params_zhdanov_general_t use params_zhdanov_species_m, only : params_zhdanov_species_t use mms_zhdanov_m, only : mms_diagnostics_zhdanov, MmsSolDens use precision_grillix_m, only : GP, MPI_GP ! From PARALLAX use screen_io_m, only : get_stdout use netcdf implicit none public :: model_zhdanov contains subroutine model_zhdanov(comm_handler) !! Runs Zhdanov model type(comm_handler_t), intent(in) :: comm_handler !! MPI communicators real (GP) :: tau !! Time integer :: isnaps !! Snapshot frame type(variable_t) :: dens_species !! Species density integer :: tstep_count !! Timestep counter type(checkpoint_monitor_t) :: checkpoint_monitor !! Main checkpoint monitor character(len=5) :: species_c !! Species counter in character format character(len=5) :: plane_c !! Plane counter in character format character(len=256) :: paramfile_zhdanov !! Parameterfile, where generic parameters for Zhdanov model !! are read from (common to all species) character(len=256) :: paramfile_zhdanov_species !! Parameterfile, where species dependend parameters for Zhdanov model !! are read from type(params_zhdanov_general_t) :: params_general !! General parameters related with ZHDANOV model type(params_zhdanov_species_t) :: params_species !! Species related parameters related with ZHDANOV model integer :: nvars type(char_arr_t), allocatable, dimension(:) :: varnames type(snapshot_t) :: snaps_zhdanov !! Snapshot handler type(karniadakis_t) :: tstep_continuity !! Time integrator for continuity equation if (is_rank_info_writer) then write(get_stdout(),*)'' write(get_stdout(),*)'---------------------------------------------------' write(get_stdout(),*)'START RUNNING ZHDANOV MODEL' write(get_stdout(),*)'---------------------------------------------------' write(get_stdout(),*)'' endif write(species_c,'(I5.5)')comm_handler%get_species() write(plane_c,'(I5.5)')comm_handler%get_plane() paramfile_zhdanov = 'params_zhdanov.in' call params_general%read_zh_controls(paramfile_zhdanov) call params_general%read_zh(paramfile_zhdanov) call params_general%read_params_timestep(paramfile_zhdanov) call params_general%display_zh(comm_handler) if (params_general%mms_on) then paramfile_zhdanov_species = '../params_zhdanov_species_'//species_c//'.in' else paramfile_zhdanov_species = 'params_zhdanov_species_'//species_c//'.in' endif call params_species%read_zh(paramfile_zhdanov_species) call params_species%display_zh(comm_handler) ! Set snapsfiles ------------------------------------------------------ if (params_general%mms_on) then nvars = 2 allocate(varnames(nvars)) varnames(2)%str = 'mms_dens' else nvars = 1 allocate(varnames(nvars)) endif varnames(1)%str = 'dens' call snaps_zhdanov%init(comm_handler, mesh, mesh, 'snapsdir_zhdanov', nvars, & varnames) ! Build initial state ------------------------------------------------- call build_zhdanov_initial_state(comm_handler, equi, mesh, map, & penalisation, & params_general, & params_species, & paramfile_zhdanov_species, & tau, isnaps, & dens_species) ! Write snapshot call write_zhdanov() ! Initialise checkpoint monitor --------------------------------------- call checkpoint_monitor%init(major_interval=params_general%tau_snaps, & major_tally_offset=isnaps) ! Initialise time-step integrators ------------------------------------ call initialise_timesteppers(comm_handler, equi, mesh, params_general,& tstep_continuity) ! Time-stepping loop -------------------------------------------------- if (is_rank_info_writer) then write(get_stdout(),*) write(get_stdout(),*)'start timestepping .....' write(get_stdout(),*) endif tstep_count = 0 do while (tau <= params_general%tau_fin + params_general%dtau) tstep_count = tstep_count + 1 call timestep_zhdanov(comm_handler, equi, equi_storage, mesh, map,& penalisation, polars, & params_general, & params_species, & tau, & dens_species, & tstep_continuity) ! Time-advancement tau=tau + params_general%dtau ! Advance checkpoint monitor call checkpoint_monitor%drive(tau, params_general%dtau) if (checkpoint_monitor%on_major_checkpoint()) then ! Update isnaps isnaps = checkpoint_monitor%get_major_tally() ! Write snapshot to snapsfile if (is_rank_info_writer) then write(get_stdout(),*)'Writing snapshot at -----------------------------------' write(get_stdout(),733)tau, tstep_count, isnaps 733 FORMAT(3X,'tau = ',ES14.6E3 /, & 3X,'tstep_count = ',I8 /, & 3X,'isnaps = ',I8 ) write(get_stdout(),*)'-------------------------------------------------------' endif ! MMS diagnostics if (params_general%mms_on) then call mms_diagnostics_zhdanov(comm_handler, equi, mesh,& tau, dens_species) endif ! Write snapshot call write_zhdanov() endif enddo if (is_rank_info_writer) then write(get_stdout(),*)'' write(get_stdout(),*)'-------------------------------------------------------' write(get_stdout(),*)'ZHDANOV MODEL RUN SUCCESFULLY' write(get_stdout(),*)'-------------------------------------------------------' write(get_stdout(),*)'' endif contains subroutine write_zhdanov() !! Writes variables of ZHDANOV model to netcdf file integer :: l real(GP) :: phi, phi_stag, x, y type(variable_t) :: mms_dens if (params_general%mms_on) then ! Compute analytic MMS solution phi = map%get_phi_cano() phi_stag = map%get_phi_stag() call mms_dens%init(comm_handler, mesh%get_n_points(), staggered = .false.) !$OMP PARALLEL PRIVATE(l, x, y) !$OMP DO do l = 1, mesh%get_n_points() x = mesh%get_x(l) y = mesh%get_y(l) mms_dens%vals(l) = MmsSolDens(equi, & comm_handler%get_species(), & x, y, phi, tau) enddo !$OMP END DO !$OMP END PARALLEL call snaps_zhdanov%write_snaps(comm_handler, isnaps, tau, & [dens_species, mms_dens]) else call snaps_zhdanov%write_snaps(comm_handler, isnaps, tau, & [dens_species]) endif end subroutine end subroutine end module