submodule (static_data_m) static_data_parbnd_immersed_s !! Handles static data related with bnd_immersed implicit none contains module subroutine static_parbnd_immersed_init(comm_handler, par_filepath, nf90_filepath) type(comm_handler_t), intent(in) :: comm_handler character(len=*), intent(in), optional :: par_filepath character(len=*), intent(in), optional :: nf90_filepath real(GP) :: dphi character(len=:), allocatable :: fpath logical :: write_to_netcdf = .false. integer :: idbg, nf90_stat, nf90_id, nf90_grpid character(len=3) :: plane_c character(len=:), allocatable :: nf90_fpath if (is_rank_info_writer) then write(get_stdout(),*)'' write(get_stdout(),*)repeat('-',80) write(get_stdout(),*)'Initializing static data: parbnd_immersed' endif if (.not.equi_is_initialized) then call handle_error('equi must be initialized at parbnd_immersed initialization', & GRILLIX_ERR_STATIC_DATA, __LINE__, __FILE__) endif if (.not.multigrid_is_initialized) then call handle_error('multigrid must be initialized at parbnd_immersed initialization', & GRILLIX_ERR_STATIC_DATA, __LINE__, __FILE__) endif if (.not.map_is_initialized) then call handle_error('map must be initialized at parbnd_immersed initialization', & GRILLIX_ERR_STATIC_DATA, __LINE__, __FILE__) endif ! Set debug output level if (enable_debug_output) then idbg = 1 else idbg = 0 endif if (present(par_filepath)) then fpath = par_filepath else fpath = static_data_parpath endif dphi = TWO_PI / comm_handler%get_nplanes() ! Select NetCDF file path and enable writing when required write(plane_c,'(I3.3)')comm_handler%get_plane() if (present(nf90_filepath)) then nf90_fpath = nf90_filepath else if (equi%is_axisymmetric()) then nf90_fpath = trim(static_data_dirpath)//'/penalisation.nc' if ((cntrl_build_parbnd_immersed == BUILD_AND_WRITE) .and. (is_rank_info_writer)) then write_to_netcdf = .true. endif else nf90_fpath = trim(static_data_dirpath)//'/penalisation_plane'//plane_c//'.nc' if (cntrl_build_parbnd_immersed == BUILD_AND_WRITE) then write_to_netcdf = .true. endif endif endif ! Associate staggered type if (equi%is_axisymmetric()) then parbnd_immersed_stag => parbnd_immersed_cano else allocate(parbnd_immersed_stag_mem) parbnd_immersed_stag => parbnd_immersed_stag_mem endif select case(cntrl_build_parbnd_immersed) case(BUILD_NOWRITE, BUILD_AND_WRITE) if (is_rank_info_writer) then write(get_stdout(),*)'Building parbnd_immersed from scratch' endif call parbnd_immersed_cano%set_parameters(dphi, filename=fpath) call parbnd_immersed_cano%build(comm_handler, equi, mesh_cano, multigrid_cano, dbgout=idbg) if (.not.equi%is_axisymmetric()) then call parbnd_immersed_stag%set_parameters(dphi, filename=fpath) call parbnd_immersed_stag%build(comm_handler, equi, mesh_stag, multigrid_stag, dbgout=idbg) endif parbnd_immersed_is_initialized = .true. case(READ) if (is_rank_info_writer) then write(get_stdout(),*)'Reading parbnd_immersed from file: ' // nf90_fpath endif nf90_stat = nf90_open(nf90_fpath, NF90_NETCDF4, nf90_id) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) if (equi%is_axisymmetric()) then call parbnd_immersed_cano%read_netcdf(nf90_id) else nf90_stat = nf90_inq_grp_ncid(nf90_id, 'canonical', nf90_grpid) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) call parbnd_immersed_cano%read_netcdf(nf90_grpid) nf90_stat = nf90_inq_grp_ncid(nf90_id, 'staggered', nf90_grpid) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) call parbnd_immersed_stag%read_netcdf(nf90_grpid) endif nf90_stat = nf90_close(nf90_id) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) parbnd_immersed_is_initialized = .true. case(NONE) if (is_rank_info_writer) then write(get_stdout(),*) 'Parbnd_immersed not built because it was not requested by the user.' endif case default call handle_error('cntrl_build_parbnd_immersed not valid', & GRILLIX_ERR_STATIC_DATA, __LINE__, __FILE__) end select if (enable_debug_output .and. parbnd_immersed_is_initialized) then call parbnd_immersed_cano%display() endif if (write_to_netcdf) then if (is_rank_info_writer) then write(get_stdout(),*)'Writing parbnd_immersed data to file (series): ' // nf90_fpath endif nf90_stat = nf90_create(nf90_fpath, NF90_NETCDF4, nf90_id) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) if (equi%is_axisymmetric()) then call parbnd_immersed_cano%write_netcdf(nf90_id) else nf90_stat = nf90_def_grp(nf90_id, 'canonical', nf90_grpid) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) call parbnd_immersed_cano%write_netcdf(nf90_grpid) nf90_stat = nf90_def_grp(nf90_id, 'staggered', nf90_grpid) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) call parbnd_immersed_stag%write_netcdf(nf90_grpid) endif nf90_stat = nf90_close(nf90_id) call handle_error_netcdf(nf90_stat, __LINE__, __FILE__) endif if (is_rank_info_writer) then write(get_stdout(),*)repeat('-',80) write(get_stdout(),*)'' endif end subroutine end submodule