diff --git a/.claude/rules/common-pitfalls.md b/.claude/rules/common-pitfalls.md index a9b69accf4..2f790a8e90 100644 --- a/.claude/rules/common-pitfalls.md +++ b/.claude/rules/common-pitfalls.md @@ -43,16 +43,27 @@ covered in `docs/documentation/contributing.md`. `toolchain/mfc/params/definitions.py`; `case_validator.py` only if physics-constrained (with a `PHYSICS_DOCS` entry). Fortran declarations and namelist bindings are auto-generated at CMake configure time — re-run cmake (or `./mfc.sh build`) after editing. -- Still manual: array variables and derived-type members (declare in - `src/*/m_global_parameters.fpp` / the relevant type), and case-optimization parameters - (`CASE_OPT_PARAMS` + the `#:else` block in `src/simulation/m_global_parameters.fpp`). - Gotcha: under `--case-optimization` those are baked into the binary and dropped from - the namelist, so changing one needs a *rebuild*, not a case edit. +- Still manual: derived-type `TYPE` member definitions in `src/common/m_derived_types.fpp`; + default-value assignments in `s_assign_default_values_to_user_inputs`; the + `CASE_OPT_EXTRA_LINES` literal in `toolchain/mfc/params/generators/fortran_gen.py` (covers `num_dims`, + `num_vels`, `weno_polyn`, `muscl_polyn`, `weno_num_stencils`, `wenojs`); + multi-variable declaration lines (`bc_x/y/z`, `x/y/z_domain`, `x/y/z_output`, post's + `G`); and MPI broadcast lists in `m_mpi_proxy`. Everything else — scalar declarations, + plain arrays (`FORTRAN_ARRAY_DIMS` table in `definitions.py`), derived-type namelist + declarations including `GPU_DECLARE` lines and Doxygen descs (`TYPED_DECLS` table in + `definitions.py`), and the simulation case-optimization declaration block — is + auto-generated at CMake configure time. + Gotcha: after editing a generator or table, force regen via `cmake_gen.py` into + `build/staging/*/` or simply reconfigure (`./mfc.sh build`) — cached builds compile + stale includes. Under `--case-optimization` the baked-in constants are dropped from the + namelist, so changing one needs a *rebuild*, not a case edit. - Runtime checks (`@:PROHIBIT`) go where they run: shared → `src/common/m_checker_common.fpp`; simulation-only → `src/simulation/m_checker.fpp`; pre/post-only → `src/{pre,post}_process/m_checker.fpp` (their `s_check_inputs` are currently empty — that IS the right place, not m_checker_common). -- Analytic ICs are compiled into the binary (syntax error = build failure, not runtime). +- Analytic ICs are compiled into the binary. Expressions are AST-validated at case load + (syntax errors and unknown variables are immediate, named errors; bare `e` is not a + variable — write `exp(1.0)`). Each IC variable maps to an `eqn_idx%…` expression in `QPVF_IDX_VARS` (`toolchain/mfc/case.py`); a new patch-settable conserved variable means updating that map AND the Fortran `eqn_idx` builder to agree — a mismatch is a silent wrong index. diff --git a/CMakeLists.txt b/CMakeLists.txt index 83bbb8fe0e..8a9c15907f 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -418,9 +418,9 @@ macro(HANDLE_SOURCES target useCommon) # Gather: # * src/[,common]/include/*.fpp - # * (if any) /include//*.fpp - file(GLOB ${target}_incs CONFIGURE_DEPENDS "${${target}_DIR}/include/*.fpp" - "${CMAKE_BINARY_DIR}/include/${target}/*.fpp") + # * generated includes from build/include// (explicit list, not a GLOB) + file(GLOB ${target}_incs CONFIGURE_DEPENDS "${${target}_DIR}/include/*.fpp") + list(APPEND ${target}_incs ${_mfc_gen_files_${target}}) if (${useCommon}) file(GLOB common_incs CONFIGURE_DEPENDS "${common_DIR}/include/*.fpp") @@ -494,6 +494,47 @@ if(_mfc_needs_regen) file(TOUCH "${_mfc_gen_stamp}") endif() +# Enumerate the 10 generated .fpp files explicitly so ninja can track them as +# build-time outputs and so HANDLE_SOURCES does not need a configure-time GLOB +# of ${CMAKE_BINARY_DIR}/include// (which fails when the dir is empty). +set(_mfc_gen_inc "${CMAKE_BINARY_DIR}/include") +set(_mfc_gen_files_pre_process + "${_mfc_gen_inc}/pre_process/generated_namelist.fpp" + "${_mfc_gen_inc}/pre_process/generated_decls.fpp" + "${_mfc_gen_inc}/pre_process/generated_constants.fpp" +) +set(_mfc_gen_files_simulation + "${_mfc_gen_inc}/simulation/generated_namelist.fpp" + "${_mfc_gen_inc}/simulation/generated_decls.fpp" + "${_mfc_gen_inc}/simulation/generated_constants.fpp" + "${_mfc_gen_inc}/simulation/generated_case_opt_decls.fpp" +) +set(_mfc_gen_files_post_process + "${_mfc_gen_inc}/post_process/generated_namelist.fpp" + "${_mfc_gen_inc}/post_process/generated_decls.fpp" + "${_mfc_gen_inc}/post_process/generated_constants.fpp" +) +set(_mfc_gen_files_syscheck) +set(_mfc_all_gen_files + ${_mfc_gen_files_pre_process} + ${_mfc_gen_files_simulation} + ${_mfc_gen_files_post_process} +) + +# Build-time regeneration: ninja re-runs cmake_gen.py and re-preprocesses +# any .fpp that includes a generated file whenever a params .py changes. +add_custom_command( + OUTPUT ${_mfc_all_gen_files} + COMMAND "${Python3_EXECUTABLE}" + "${CMAKE_CURRENT_SOURCE_DIR}/toolchain/mfc/params/generators/cmake_gen.py" + "${CMAKE_BINARY_DIR}" + DEPENDS ${_mfc_gen_inputs} + WORKING_DIRECTORY "${CMAKE_CURRENT_SOURCE_DIR}" + COMMENT "Regenerating Fortran parameter includes" + VERBATIM +) +add_custom_target(mfc_params_gen DEPENDS ${_mfc_all_gen_files}) + HANDLE_SOURCES(pre_process ON) HANDLE_SOURCES(simulation ON) HANDLE_SOURCES(post_process ON) diff --git a/docs/documentation/contributing.md b/docs/documentation/contributing.md index 6bec8f738d..a79c52006c 100644 --- a/docs/documentation/contributing.md +++ b/docs/documentation/contributing.md @@ -276,35 +276,27 @@ def check_my_feature(self): If your check enforces a physics constraint, also add a `PHYSICS_DOCS` entry (see [How to Document Physics Constraints](#how-to-document-physics-constraints) below). -**Step 5: Declare in Fortran** (`src//m_global_parameters.fpp`) +**Step 5: Fortran declaration and namelist binding (auto-generated)** -Add the variable declaration in the appropriate target's global parameters module. Choose the target(s) where the parameter is used: +Scalar declarations, GPU declare lines, Doxygen descriptions, and namelist bindings are +auto-generated at CMake configure time from the `TYPED_DECLS` and `FORTRAN_ARRAY_DIMS` +tables in `toolchain/mfc/params/definitions.py`. For a plain scalar registered with +`_r()` / `_nv()` above, no manual Fortran edit is needed — reconfigure (`./mfc.sh build`) +and the generated include in each target's `m_global_parameters.fpp` is updated +automatically. -- `src/pre_process/m_global_parameters.fpp` -- `src/simulation/m_global_parameters.fpp` -- `src/post_process/m_global_parameters.fpp` +Still manual (not auto-generated): -```fortran -real(wp) :: my_param !< Description of the parameter -``` - -If the parameter is used in GPU kernels, add a GPU declaration: - -```fortran -$:GPU_DECLARE(create='[my_param]') -``` - -**Step 6: Add to Fortran namelist** (`src//m_start_up.fpp`) - -Add the parameter name to the `namelist /user_inputs/` declaration: - -```fortran -namelist /user_inputs/ ... , my_param, ... -``` +- `TYPE` member definitions inside derived types in `src/common/m_derived_types.fpp` +- Default-value assignments in `s_assign_default_values_to_user_inputs` +- Multi-variable declaration lines (`bc_x/y/z`, `x/y/z_domain`, `x/y/z_output`) +- MPI broadcast lists in `src/*/m_mpi_proxy.fpp` +- `CASE_OPT_EXTRA_LINES` in `toolchain/mfc/params/generators/fortran_gen.py` for case-optimization constants -The toolchain writes the parameter to the input file and Fortran reads it via this namelist. No other I/O code is needed. +After editing any generator or table, force regen by reconfiguring (`./mfc.sh build`) — +cached builds compile stale includes. -**Step 7: Use in Fortran code** +**Step 6: Use in Fortran code** Reference `my_param` anywhere in the target's modules. It is available as a global after the namelist is read at startup. diff --git a/src/common/m_constants.fpp b/src/common/m_constants.fpp index 5452c13a21..18d4553c1c 100644 --- a/src/common/m_constants.fpp +++ b/src/common/m_constants.fpp @@ -50,10 +50,7 @@ module m_constants real(wp), parameter :: broadband_spectral_level_constant = 20._wp !> The spectral level constant to correct the magnitude at each frequency to ensure the source is overall broadband real(wp), parameter :: broadband_spectral_level_growth_rate = 10._wp - ! Reconstruction Types - integer, parameter :: WENO_TYPE = 1 !< Using WENO for reconstruction type - integer, parameter :: MUSCL_TYPE = 2 !< Using MUSCL for reconstruction type - ! Interface Compression + ! Reconstruction Types Interface Compression real(wp), parameter :: dflt_ic_eps = 1e-4_wp !< Ensure compression is only applied to surface cells in THINC real(wp), parameter :: dflt_ic_beta = 1.6_wp !< Sharpness parameter's default value used in THINC real(wp), parameter :: moncon_cutoff = 1e-8_wp !< Monotonicity constraint's limiter to prevent extremas in THINC diff --git a/src/common/m_helper_basic.fpp b/src/common/m_helper_basic.fpp index ccd46145e1..be1e89c398 100644 --- a/src/common/m_helper_basic.fpp +++ b/src/common/m_helper_basic.fpp @@ -9,6 +9,7 @@ module m_helper_basic use m_derived_types use m_precision_select + use m_constants, only: recon_type_weno, recon_type_muscl implicit none @@ -125,13 +126,13 @@ contains if (igr) then buff_size = (igr_order - 1)/2 + 2 - else if (recon_type == WENO_TYPE) then + else if (recon_type == recon_type_weno) then if (viscous) then buff_size = 2*weno_polyn + 2 else buff_size = weno_polyn + 2 end if - else if (recon_type == MUSCL_TYPE) then + else if (recon_type == recon_type_muscl) then buff_size = muscl_polyn + 2 end if diff --git a/src/common/m_mpi_common.fpp b/src/common/m_mpi_common.fpp index 9cfdd82e14..b509e1f6ea 100644 --- a/src/common/m_mpi_common.fpp +++ b/src/common/m_mpi_common.fpp @@ -17,6 +17,7 @@ module m_mpi_common use m_helper use ieee_arithmetic use m_nvtx + use m_constants, only: recon_type_weno, format_silo implicit none @@ -1031,7 +1032,7 @@ contains integer :: i, j !< Generic loop iterators integer :: ierr !< Generic flag used to identify and report MPI errors - if (recon_type == WENO_TYPE) then + if (recon_type == recon_type_weno) then recon_order = weno_order else recon_order = muscl_order @@ -1194,14 +1195,14 @@ contains #ifdef MFC_POST_PROCESS ! Ghost zone at the beginning - if (proc_coords(3) > 0 .and. format == 1) then + if (proc_coords(3) > 0 .and. format == format_silo) then offset_z%beg = 2 else offset_z%beg = 0 end if ! Ghost zone at the end - if (proc_coords(3) < num_procs_z - 1 .and. format == 1) then + if (proc_coords(3) < num_procs_z - 1 .and. format == format_silo) then offset_z%end = 2 else offset_z%end = 0 @@ -1306,14 +1307,14 @@ contains #ifdef MFC_POST_PROCESS ! Ghost zone at the beginning - if (proc_coords(2) > 0 .and. format == 1) then + if (proc_coords(2) > 0 .and. format == format_silo) then offset_y%beg = 2 else offset_y%beg = 0 end if ! Ghost zone at the end - if (proc_coords(2) < num_procs_y - 1 .and. format == 1) then + if (proc_coords(2) < num_procs_y - 1 .and. format == format_silo) then offset_y%end = 2 else offset_y%end = 0 @@ -1389,14 +1390,14 @@ contains #ifdef MFC_POST_PROCESS ! Ghost zone at the beginning - if (proc_coords(1) > 0 .and. format == 1) then + if (proc_coords(1) > 0 .and. format == format_silo) then offset_x%beg = 2 else offset_x%beg = 0 end if ! Ghost zone at the end - if (proc_coords(1) < num_procs_x - 1 .and. format == 1) then + if (proc_coords(1) < num_procs_x - 1 .and. format == format_silo) then offset_x%end = 2 else offset_x%end = 0 diff --git a/src/common/m_variables_conversion.fpp b/src/common/m_variables_conversion.fpp index 97c9dd0c4a..db7f80faf2 100644 --- a/src/common/m_variables_conversion.fpp +++ b/src/common/m_variables_conversion.fpp @@ -14,7 +14,7 @@ module m_variables_conversion use m_helper_basic use m_helper use m_constants, only: riemann_solver_hll, riemann_solver_hlld, model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, & - & model_eqns_4eq + & model_eqns_4eq, avg_state_roe use m_thermochem, only: num_species, get_temperature, get_pressure, gas_constant, get_mixture_molecular_weight, & & get_mixture_energy_mass @@ -1267,7 +1267,7 @@ contains integer :: q if (chemistry) then ! Reacting mixture sound speed - if (avg_state == 1 .and. abs(c_c) > verysmall) then + if (avg_state == avg_state_roe .and. abs(c_c) > verysmall) then c = sqrt(c_c - (gamma - 1.0_wp)*(vel_sum - H)) else c = sqrt((1.0_wp + 1.0_wp/gamma)*pres/rho) diff --git a/src/post_process/m_data_output.fpp b/src/post_process/m_data_output.fpp index 1a189408ac..3b40026a7d 100644 --- a/src/post_process/m_data_output.fpp +++ b/src/post_process/m_data_output.fpp @@ -12,7 +12,7 @@ module m_data_output use m_compile_specific use m_helper use m_variables_conversion - use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq + use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, format_silo, format_binary, precision_single implicit none @@ -87,7 +87,7 @@ contains allocate (cyl_q_sf(-offset_y%beg:n + offset_y%end,-offset_z%beg:p + offset_z%end,-offset_x%beg:m + offset_x%end)) end if - if (precision == 1) then + if (precision == precision_single) then allocate (q_sf_s(-offset_x%beg:m + offset_x%end,-offset_y%beg:n + offset_y%end,-offset_z%beg:p + offset_z%end)) if (grid_geometry == 3) then allocate (cyl_q_sf_s(-offset_y%beg:n + offset_y%end,-offset_z%beg:p + offset_z%end,-offset_x%beg:m + offset_x%end)) @@ -96,7 +96,7 @@ contains if (n == 0) then allocate (q_root_sf(0:m_root,0:0,0:0)) - if (precision == 1) then + if (precision == precision_single) then allocate (q_root_sf_s(0:m_root,0:0,0:0)) end if end if @@ -104,7 +104,7 @@ contains ! Allocating the spatial and data extents and also the variables for the offsets and the one bookkeeping the number of ! cell-boundaries in each active coordinate direction. Note that all these variables are only needed by the Silo-HDF5 format ! for multidimensional data. - if (format == 1) then + if (format == format_silo) then allocate (data_extents(1:2,0:num_procs - 1)) if (p > 0) then @@ -129,7 +129,7 @@ contains ! results are now transferred to the local variables of this module when they are required by the Silo-HDF5 format, for ! multidimensional data sets. With the same, latter, requirements, the variables bookkeeping the number of cell-boundaries ! in each active coordinate direction are also set here. - if (format == 1) then + if (format == format_silo) then if (p > 0) then if (grid_geometry == 3) then lo_offset(:) = (/offset_y%beg, offset_z%beg, offset_x%beg/) @@ -158,7 +158,7 @@ contains end if end if - if (format == 1) then + if (format == format_silo) then dbdir = trim(case_dir) // '/silo_hdf5' write (proc_rank_dir, '(A,I0)') '/p', proc_rank @@ -225,12 +225,12 @@ contains ! Contrary to the Silo-HDF5 database format, handles of the Binary database master/root and slave/local process files are ! perfectly static throughout post-process. Hence, they are set here so that they do not have to be repetitively computed in ! later procedures. - if (format == 2) then + if (format == format_binary) then if (n == 0 .and. proc_rank == 0) dbroot = 2 dbfile = 1 end if - if (format == 2) then + if (format == format_binary) then dbvars = 0 if ((model_eqns == model_eqns_5eq) .or. (model_eqns == model_eqns_6eq)) then @@ -357,7 +357,7 @@ contains character(LEN=len_trim(case_dir) + 3*name_len) :: file_loc integer :: ierr - if (format == 1) then + if (format == format_silo) then write (file_loc, '(A,I0,A)') '/', t_step, '.silo' file_loc = trim(proc_rank_dir) // trim(file_loc) @@ -450,7 +450,7 @@ contains integer :: i integer :: ierr - if (format == 1) then + if (format == format_silo) then ! For multidimensional data sets, the spatial extents of all of the grid(s) handled by the local processor(s) are ! recorded so that they may be written, by root processor, to the formatted database master file. if (num_procs > 1) then @@ -516,11 +516,11 @@ contains & DB_DOUBLE, DB_COLLINEAR, optlist, ierr) err = DBFREEOPTLIST(optlist) end if - else if (format == 2) then + else if (format == format_binary) then ! Multidimensional local grid data is written to the formatted database slave file. Recall that no master file to ! maintained in multidimensions. if (p > 0) then - if (precision == 1) then + if (precision == precision_single) then write (dbfile) real(x_cb, sp), real(y_cb, sp), real(z_cb, sp) else if (output_partial_domain) then @@ -531,7 +531,7 @@ contains end if end if else if (n > 0) then - if (precision == 1) then + if (precision == precision_single) then write (dbfile) real(x_cb, sp), real(y_cb, sp) else if (output_partial_domain) then @@ -544,7 +544,7 @@ contains ! One-dimensional local grid data is written to the formatted database slave file. In addition, the local grid data ! is put together by the root process and written to the master file. else - if (precision == 1) then + if (precision == precision_single) then write (dbfile) real(x_cb, sp) else if (output_partial_domain) then @@ -561,7 +561,7 @@ contains end if if (proc_rank == 0) then - if (precision == 1) then + if (precision == precision_single) then write (dbroot) real(x_root_cb, wp) else if (output_partial_domain) then @@ -588,7 +588,7 @@ contains integer :: i, j, k integer :: ierr - if (format == 1) then + if (format == format_silo) then ! Determining the extents of the flow variable on each local process and gathering all this information on root process if (num_procs > 1) then call s_mpi_gather_data_extents(q_sf, data_extents) @@ -613,7 +613,7 @@ contains end if if (wp == dp) then - if (precision == 1) then + if (precision == precision_single) then do i = -offset_x%beg, m + offset_x%end do j = -offset_y%beg, n + offset_y%end do k = -offset_z%beg, p + offset_z%end @@ -682,7 +682,7 @@ contains else ! Writing the name of the flow variable and its data, associated with the local processor, to the formatted database ! slave file - if (precision == 1) then + if (precision == precision_single) then write (dbfile) varname, real(q_sf, wp) else write (dbfile) varname, q_sf @@ -698,7 +698,7 @@ contains end if if (proc_rank == 0) then - if (precision == 1) then + if (precision == precision_single) then write (dbroot) varname, real(q_root_sf, wp) else write (dbroot) varname, q_root_sf @@ -1496,7 +1496,7 @@ contains integer :: ierr - if (format == 1) then + if (format == format_silo) then ierr = DBCLOSE(dbfile) if (proc_rank == 0) ierr = DBCLOSE(dbroot) else @@ -1532,7 +1532,7 @@ contains ! Deallocating spatial and data extents and also the variables for the offsets and the one bookkeeping the number of ! cell-boundaries in each active coordinate direction. Note that all these variables were only needed by Silo-HDF5 format ! for multidimensional data. - if (format == 1) then + if (format == format_silo) then deallocate (spatial_extents) deallocate (data_extents) deallocate (lo_offset) diff --git a/src/post_process/m_global_parameters.fpp b/src/post_process/m_global_parameters.fpp index 305178ec0f..221042645e 100644 --- a/src/post_process/m_global_parameters.fpp +++ b/src/post_process/m_global_parameters.fpp @@ -14,7 +14,8 @@ module m_global_parameters use m_derived_types use m_helper_basic use m_thermochem, only: num_species, species_names - use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, model_eqns_4eq + use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, model_eqns_4eq, format_silo, precision_single, & + & recon_type_weno implicit none @@ -121,10 +122,8 @@ module m_global_parameters integer :: mpi_info_int !> @} - type(physical_parameters), dimension(num_fluids_max) :: fluid_pp !< Stiffened gas EOS parameters and Reynolds numbers per fluid - ! Subgrid Bubble Parameters - type(subgrid_bubble_physical_parameters) :: bub_pp - real(wp), allocatable, dimension(:) :: adv !< Advection variables + ! fluid_pp, bub_pp: auto-generated in generated_decls.fpp + real(wp), allocatable, dimension(:) :: adv !< Advection variables ! Formatted Database File(s) Structure Parameters type(bounds_info) :: x_output, y_output, z_output !< Portion of domain to output for post-processing @@ -188,7 +187,7 @@ contains ! Simulation algorithm parameters model_eqns = dflt_int num_fluids = dflt_int - recon_type = WENO_TYPE + recon_type = recon_type_weno weno_order = dflt_int muscl_order = dflt_int mixture_err = .false. @@ -668,7 +667,7 @@ contains ! Size of the ghost zone layer is non-zero only when post-processing the raw simulation data of a parallel multidimensional ! computation in the Silo-HDF5 format. If this is the case, one must also verify whether the raw simulation data is 2D or ! 3D. In the 2D case, size of the z-coordinate direction ghost zone layer must be zeroed out. - if (num_procs == 1 .or. format /= 1) then + if (num_procs == 1 .or. format /= format_silo) then offset_x%beg = 0 offset_x%end = 0 offset_y%beg = 0 @@ -733,7 +732,7 @@ contains allocate (x_root_cb(-1:m_root)) allocate (x_root_cc(0:m_root)) - if (precision == 1) then + if (precision == precision_single) then allocate (x_root_cc_s(0:m_root)) end if end if diff --git a/src/post_process/m_mpi_proxy.fpp b/src/post_process/m_mpi_proxy.fpp index e6e3769ee1..1636539d7a 100644 --- a/src/post_process/m_mpi_proxy.fpp +++ b/src/post_process/m_mpi_proxy.fpp @@ -13,6 +13,7 @@ module m_mpi_proxy use m_global_parameters use m_mpi_common use ieee_arithmetic + use m_constants, only: format_silo implicit none @@ -35,7 +36,7 @@ contains ! procedures. Note that these are only needed for either multidimensional runs that utilize the Silo database file format or ! for 1D simulations. - if ((format == 1 .and. n > 0) .or. n == 0) then + if ((format == format_silo .and. n > 0) .or. n == 0) then allocate (recvcounts(0:num_procs - 1)) allocate (displs(0:num_procs - 1)) @@ -249,7 +250,7 @@ contains integer :: ierr !< Generic flag used to identify and report MPI errors ! Silo-HDF5 database format - if (format == 1) then + if (format == format_silo) then call MPI_GATHERV(x_cc(0), m + 1, mpi_p, x_root_cc(0), recvcounts, displs, mpi_p, 0, MPI_COMM_WORLD, ierr) ! Binary database format @@ -319,7 +320,7 @@ contains #ifdef MFC_MPI ! Deallocating the receive counts and the displacement vector variables used in variable-gather communication procedures - if ((format == 1 .and. n > 0) .or. n == 0) then + if ((format == format_silo .and. n > 0) .or. n == 0) then deallocate (recvcounts) deallocate (displs) end if diff --git a/src/pre_process/m_data_output.fpp b/src/pre_process/m_data_output.fpp index 9d58377214..c0dd4c906a 100644 --- a/src/pre_process/m_data_output.fpp +++ b/src/pre_process/m_data_output.fpp @@ -21,7 +21,7 @@ module m_data_output use m_boundary_conditions use m_thermochem, only: species_names use m_helper - use m_constants, only: model_eqns_5eq + use m_constants, only: model_eqns_5eq, precision_single implicit none @@ -147,7 +147,7 @@ contains pi_inf = pi_infs(1) qv = qvs(1) - if (precision == 1) then + if (precision == precision_single) then FMT = "(2F30.3)" else FMT = "(2F40.14)" @@ -272,7 +272,7 @@ contains end if end if - if (precision == 1) then + if (precision == precision_single) then FMT = "(3F30.7)" else FMT = "(3F40.14)" @@ -323,7 +323,7 @@ contains end if end if - if (precision == 1) then + if (precision == precision_single) then FMT = "(4F30.7)" else FMT = "(4F40.14)" diff --git a/src/pre_process/m_global_parameters.fpp b/src/pre_process/m_global_parameters.fpp index 3c48bdf3e9..a5829e5348 100644 --- a/src/pre_process/m_global_parameters.fpp +++ b/src/pre_process/m_global_parameters.fpp @@ -14,7 +14,7 @@ module m_global_parameters use m_derived_types ! Definitions of the derived types use m_helper_basic ! Functions to compare floating point numbers use m_thermochem, only: num_species - use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, model_eqns_4eq + use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, model_eqns_4eq, recon_type_weno, recon_type_muscl implicit none @@ -58,15 +58,15 @@ module m_global_parameters type(int_bounds_info) :: idwint(1:3) ! Cell indices (InDices With BUFFer): includes buffer except in pre_process - type(int_bounds_info) :: idwbuff(1:3) - type(int_bounds_info) :: bc_x, bc_y, bc_z !< Boundary conditions in the x-, y- and z-coordinate directions - integer :: shear_num !< Number of shear stress components - integer, dimension(3) :: shear_indices !< Indices of the stress components that represent shear stress - integer :: shear_BC_flip_num !< Number of shear stress components to reflect for boundary conditions - integer, dimension(3, 2) :: shear_BC_flip_indices !< Shear stress BC reflection indices (1:3, 1:shear_BC_flip_num) - type(simplex_noise_params) :: simplex_params - - ! Perturb density of surrounding air so as to break symmetry of grid fluid_rho: auto-generated in generated_decls.fpp + type(int_bounds_info) :: idwbuff(1:3) + type(int_bounds_info) :: bc_x, bc_y, bc_z !< Boundary conditions in the x-, y- and z-coordinate directions + integer :: shear_num !< Number of shear stress components + integer, dimension(3) :: shear_indices !< Indices of the stress components that represent shear stress + integer :: shear_BC_flip_num !< Number of shear stress components to reflect for boundary conditions + integer, dimension(3, 2) :: shear_BC_flip_indices !< Shear stress BC reflection indices (1:3, 1:shear_BC_flip_num) + ! simplex_params: auto-generated in generated_decls.fpp + + ! fluid_rho (perturbs surrounding-air density to break grid symmetry): auto-generated in generated_decls.fpp integer, allocatable, dimension(:) :: proc_coords !< Processor coordinates in MPI_CART_COMM integer, allocatable, dimension(:) :: start_idx !< Starting cell-center index of local processor in global grid #ifdef MFC_MPI @@ -75,16 +75,11 @@ module m_global_parameters integer :: mpi_info_int !< MPI info for parallel IO with Lustre file systems #endif - ! Initial Condition Parameters - type(ic_patch_parameters), dimension(num_patches_max) :: patch_icpp !< IC patch parameters (max: num_patches_max) - type(bc_patch_parameters), dimension(num_bc_patches_max) :: patch_bc !< Boundary condition patch parameters - logical :: bc_io !< whether or not to save BC data + ! Initial Condition Parameters patch_icpp, patch_bc: auto-generated in generated_decls.fpp + logical :: bc_io !< whether or not to save BC data - ! Fluids Physical Parameters - type(physical_parameters), dimension(num_fluids_max) :: fluid_pp !< Stiffened gas EOS parameters and Reynolds numbers per fluid - ! Subgrid Bubble Parameters - type(subgrid_bubble_physical_parameters) :: bub_pp - type(chemistry_parameters) :: chem_params + ! Fluids Physical Parameters fluid_pp, bub_pp: auto-generated in generated_decls.fpp + type(chemistry_parameters) :: chem_params !> @name Bubble modeling !> @{ real(wp) :: Eu @@ -92,12 +87,7 @@ module m_global_parameters integer :: nmom !< Number of carried moments !> @} - !> @name Immersed Boundaries - !> @{ - type(ib_patch_parameters), dimension(num_ib_patches_max_namelist) :: patch_ib !< Immersed boundary patch parameters - type(ib_airfoil_parameters), dimension(num_ib_airfoils_max) :: ib_airfoil !< Per-airfoil NACA parameters - type(ib_stl_parameters), dimension(num_stl_models_max) :: stl_models !< Per-STL model parameters (namelist) - !> @} + ! patch_ib, ib_airfoil, stl_models: auto-generated in generated_decls.fpp !> @name Non-polytropic bubble gas compression !> @{ @@ -171,7 +161,7 @@ contains palpha_eps = dflt_real ptgalpha_eps = dflt_real num_fluids = dflt_int - recon_type = WENO_TYPE + recon_type = recon_type_weno weno_order = dflt_int igr = .false. igr_order = dflt_int @@ -460,9 +450,9 @@ contains integer :: i, j, fac - if (recon_type == WENO_TYPE) then + if (recon_type == recon_type_weno) then weno_polyn = (weno_order - 1)/2 - else if (recon_type == MUSCL_TYPE) then + else if (recon_type == recon_type_muscl) then muscl_polyn = muscl_order end if diff --git a/src/simulation/include/inline_riemann.fpp b/src/simulation/include/inline_riemann.fpp index 70734870ad..b5ededfaa1 100644 --- a/src/simulation/include/inline_riemann.fpp +++ b/src/simulation/include/inline_riemann.fpp @@ -59,11 +59,11 @@ #:enddef roe_avg #:def compute_average_state() - if (avg_state == 1) then + if (avg_state == avg_state_roe) then @:roe_avg() end if - if (avg_state == 2) then + if (avg_state == avg_state_arithmetic) then @:arithmetic_avg() end if #:enddef compute_average_state diff --git a/src/simulation/m_bubbles.fpp b/src/simulation/m_bubbles.fpp index 1881b128b4..2a7483735f 100644 --- a/src/simulation/m_bubbles.fpp +++ b/src/simulation/m_bubbles.fpp @@ -12,6 +12,7 @@ module m_bubbles use m_mpi_proxy use m_variables_conversion use m_helper_basic + use m_constants, only: bubble_model_gilmore, bubble_model_keller_miksis, bubble_model_rayleigh_plesset implicit none @@ -32,7 +33,7 @@ contains real(wp) :: fCpbw, fCpinf, fCpinf_dot, fH, fHdot, c_gas, c_liquid real(wp) :: f_rddot - if (bubble_model == 1) then + if (bubble_model == bubble_model_gilmore) then ! Gilmore bubbles fCpinf = fP - Eu fCpbw = f_cpbw(fR0, fR, fV, fpb) @@ -41,7 +42,7 @@ contains fCpinf_dot = f_cpinfdot(fRho, fP, alf, fntait, fBtait, f_bub_adv_src, f_divu) fHdot = f_Hdot(fCpbw, fCpinf, fCpinf_dot, fntait, fBtait, fR, fV, fR0, fpbdot) f_rddot = f_rddot_G(fCpbw, fR, fV, fH, fHdot, c_gas, fntait, fBtait) - else if (bubble_model == 2) then + else if (bubble_model == bubble_model_keller_miksis) then ! Keller-Miksis bubbles fCpinf = fP fCpbw = f_cpbw_KM(fR0, fR, fV, fpb) @@ -51,7 +52,7 @@ contains c_liquid = fCson end if f_rddot = f_rddot_KM(fpbdot, fCpinf, fCpbw, fRho, fR, fV, fR0, c_liquid) - else if (bubble_model == 3) then + else if (bubble_model == bubble_model_rayleigh_plesset) then ! Rayleigh-Plesset bubbles fCpbw = f_cpbw_KM(fR0, fR, fV, fpb) f_rddot = f_rddot_RP(fP, fRho, fR, fV, fCpbw) diff --git a/src/simulation/m_bubbles_EL.fpp b/src/simulation/m_bubbles_EL.fpp index ab6cc6c461..4a84874119 100644 --- a/src/simulation/m_bubbles_EL.fpp +++ b/src/simulation/m_bubbles_EL.fpp @@ -17,7 +17,7 @@ module m_bubbles_EL use m_helper_basic use m_sim_helpers use m_helper - use m_constants, only: time_stepper_rk1, time_stepper_rk2, time_stepper_rk3 + use m_constants, only: time_stepper_rk1, time_stepper_rk2, time_stepper_rk3, precision_single implicit none @@ -1281,7 +1281,7 @@ contains file_loc = trim(case_dir) // '/D/' // trim(file_loc) inquire (FILE=trim(file_loc), EXIST=file_exist) - if (precision == 1) then + if (precision == precision_single) then FMT = "(A16,A14,8A16)" else FMT = "(A24,A14,8A24)" @@ -1295,7 +1295,7 @@ contains open (11, FILE=trim(file_loc), form='formatted', position='append') end if - if (precision == 1) then + if (precision == precision_single) then FMT = "(F16.8,I14,8F16.8)" else FMT = "(F24.16,I14,8F24.16)" @@ -1543,7 +1543,7 @@ contains $:GPU_UPDATE(host='[Rmax_glb, Rmin_glb]') - if (precision == 1) then + if (precision == precision_single) then FMT = "(A10,A14,5A16)" else FMT = "(A10,A14,5A24)" @@ -1552,7 +1552,7 @@ contains open (13, FILE=trim(file_loc), form='formatted', position='rewind') write (13, FMT) 'proc_rank', 'particleID', 'x', 'y', 'z', 'Rmax_glb', 'Rmin_glb' - if (precision == 1) then + if (precision == precision_single) then FMT = "(I10,I14,5F16.8)" else FMT = "(I10,I14,5F24.16)" diff --git a/src/simulation/m_cbc.fpp b/src/simulation/m_cbc.fpp index 7241c4ea6b..d534169504 100644 --- a/src/simulation/m_cbc.fpp +++ b/src/simulation/m_cbc.fpp @@ -12,7 +12,7 @@ module m_cbc use m_global_parameters use m_variables_conversion use m_compute_cbc - use m_constants, only: riemann_solver_hll, model_eqns_gamma_law + use m_constants, only: riemann_solver_hll, model_eqns_gamma_law, recon_type_weno, recon_type_muscl use m_thermochem, only: get_mixture_energy_mass, get_mixture_specific_heat_cv_mass, get_mixture_specific_heat_cp_mass, & & gas_constant, get_mixture_molecular_weight, get_species_enthalpies_rt, molecular_weights, get_species_specific_heats_r, & & get_mole_fractions, get_species_specific_heats_r @@ -189,12 +189,12 @@ contains ! Allocating the cell-width distribution in the s-direction @:ALLOCATE(ds(0:buff_size)) - if (recon_type == WENO_TYPE) then + if (recon_type == recon_type_weno) then idx1%beg = 0 idx1%end = weno_polyn - 1 idx2%beg = 0 idx2%end = weno_order - 3 - else if (recon_type == MUSCL_TYPE) then + else if (recon_type == recon_type_muscl) then idx1%beg = 0 idx1%end = muscl_polyn idx2%beg = 0 @@ -350,7 +350,7 @@ contains ! Computing CBC1 Coefficients #:for CBC_DIR, XYZ in [(1, 'x'), (2, 'y'), (3, 'z')] - if (cbc_dir_in == ${CBC_DIR}$ .and. recon_type == WENO_TYPE) then + if (cbc_dir_in == ${CBC_DIR}$ .and. recon_type == recon_type_weno) then if (weno_order == 1) then fd_coef_${XYZ}$ (:,cbc_loc_in) = 0._wp fd_coef_${XYZ}$ (0, cbc_loc_in) = -2._wp/(ds(0) + ds(1)) @@ -524,7 +524,7 @@ contains call s_associate_cbc_coefficients_pointers(cbc_dir, cbc_loc) #:for CBC_DIR, XYZ in [(1, 'x'), (2, 'y'), (3, 'z')] - if (cbc_dir == ${CBC_DIR}$ .and. recon_type == WENO_TYPE) then + if (cbc_dir == ${CBC_DIR}$ .and. recon_type == recon_type_weno) then ! PI2 of flux_rs_vf and flux_src_rs_vf at j = 1/2 if (weno_order == 3) then call s_convert_primitive_to_flux_variables(q_prim_rs${XYZ}$_vf, F_rs${XYZ}$_vf, F_src_rs${XYZ}$_vf, is1, is2, & diff --git a/src/simulation/m_checker.fpp b/src/simulation/m_checker.fpp index 2ed8aca3c3..8538e273b7 100644 --- a/src/simulation/m_checker.fpp +++ b/src/simulation/m_checker.fpp @@ -12,6 +12,7 @@ module m_checker use m_mpi_proxy use m_helper use m_helper_basic + use m_constants, only: recon_type_weno, recon_type_muscl, muscl_order_first_order implicit none @@ -27,9 +28,9 @@ contains if (igr) then call s_check_inputs_nvidia_uvm else - if (recon_type == WENO_TYPE) then + if (recon_type == recon_type_weno) then call s_check_inputs_weno - else if (recon_type == MUSCL_TYPE) then + else if (recon_type == recon_type_muscl) then call s_check_inputs_muscl end if end if @@ -84,7 +85,7 @@ contains @:PROHIBIT(p + 1 < min(1, p)*num_stcls_min*muscl_order, & & "For 3D simulation, p must be greater than or equal to (num_stcls_min*muscl_order - 1), whose value is " & & // trim(numStr)) - @:PROHIBIT(muscl_order == 1 .and. int_comp > 0, & + @:PROHIBIT(muscl_order == muscl_order_first_order .and. int_comp > 0, & & "int_comp requires muscl_order >= 2 (muscl_order=1 leaves the reconstruction workspace uninitialised)") end subroutine s_check_inputs_muscl diff --git a/src/simulation/m_data_output.fpp b/src/simulation/m_data_output.fpp index 65289e3cc7..255e3a4dd3 100644 --- a/src/simulation/m_data_output.fpp +++ b/src/simulation/m_data_output.fpp @@ -19,7 +19,7 @@ module m_data_output use m_delay_file_access use m_ibm use m_boundary_common - use m_constants, only: model_eqns_5eq, model_eqns_4eq + use m_constants, only: model_eqns_5eq, model_eqns_4eq, precision_single implicit none @@ -381,7 +381,7 @@ contains pi_inf = pi_infs(1) qv = qvs(1) - if (precision == 1) then + if (precision == precision_single) then FMT = "(2F30.3)" else FMT = "(2F40.14)" @@ -462,7 +462,7 @@ contains end if end if - if (precision == 1) then + if (precision == precision_single) then FMT = "(3F30.7)" else FMT = "(3F40.14)" @@ -546,7 +546,7 @@ contains end if end if - if (precision == 1) then + if (precision == precision_single) then FMT = "(4F30.7)" else FMT = "(4F40.14)" diff --git a/src/simulation/m_global_parameters.fpp b/src/simulation/m_global_parameters.fpp index 479e62728e..feed63d798 100644 --- a/src/simulation/m_global_parameters.fpp +++ b/src/simulation/m_global_parameters.fpp @@ -14,7 +14,7 @@ module m_global_parameters use m_derived_types use m_helper_basic - use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, model_eqns_4eq + use m_constants, only: model_eqns_gamma_law, model_eqns_5eq, model_eqns_6eq, model_eqns_4eq, recon_type_weno, recon_type_muscl ! $:USE_GPU_MODULE() implicit none @@ -65,57 +65,7 @@ module m_global_parameters logical :: cfl_dt ! Simulation Algorithm Parameters - #:if MFC_CASE_OPTIMIZATION - integer, parameter :: num_dims = ${num_dims}$ !< Number of spatial dimensions - integer, parameter :: num_vels = ${num_vels}$ !< Number of velocity components (different from num_dims for mhd) - #:else - integer :: num_dims !< Number of spatial dimensions - integer :: num_vels !< Number of velocity components (different from num_dims for mhd) - #:endif - #:if MFC_CASE_OPTIMIZATION - integer, parameter :: recon_type = ${recon_type}$ !< Reconstruction type - integer, parameter :: weno_polyn = ${weno_polyn}$ !< Degree of the WENO polynomials (polyn) - integer, parameter :: muscl_polyn = ${muscl_polyn}$ !< Degree of the MUSCL polynomials (polyn) - integer, parameter :: weno_order = ${weno_order}$ !< Order of the WENO reconstruction - integer, parameter :: muscl_order = ${muscl_order}$ !< Order of the MUSCL order - !> Number of stencils for WENO reconstruction (only different from weno_polyn for TENO(>5)) - integer, parameter :: weno_num_stencils = ${weno_num_stencils}$ - integer, parameter :: muscl_lim = ${muscl_lim}$ !< MUSCL Limiter - integer, parameter :: num_fluids = ${num_fluids}$ !< number of fluids in the simulation - logical, parameter :: wenojs = (${wenojs}$ /= 0) !< WENO-JS (default) - logical, parameter :: mapped_weno = (${mapped_weno}$ /= 0) !< WENO-M (WENO with mapping of nonlinear weights) - logical, parameter :: wenoz = (${wenoz}$ /= 0) !< WENO-Z - logical, parameter :: teno = (${teno}$ /= 0) !< TENO (Targeted ENO) - real(wp), parameter :: wenoz_q = ${wenoz_q}$ !< Power constant for WENO-Z - logical, parameter :: mhd = (${mhd}$ /= 0) !< Magnetohydrodynamics - logical, parameter :: relativity = (${relativity}$ /= 0) !< Relativity (only for MHD) - integer, parameter :: igr_iter_solver = ${igr_iter_solver}$ !< IGR elliptic solver - integer, parameter :: igr_order = ${igr_order}$ !< Reconstruction order for IGR - logical, parameter :: igr = (${igr}$ /= 0) !< use information geometric regularization - logical, parameter :: igr_pres_lim = (${igr_pres_lim}$ /= 0) !< Limit to positive pressures for IGR - logical, parameter :: viscous = (${viscous}$ /= 0) !< Viscous effects - #:else - integer :: recon_type - integer :: weno_polyn - integer :: muscl_polyn - integer :: weno_order - integer :: muscl_order - integer :: weno_num_stencils - integer :: muscl_lim - integer :: num_fluids - logical :: wenojs - logical :: mapped_weno - logical :: wenoz - logical :: teno - real(wp) :: wenoz_q - logical :: mhd - logical :: relativity - integer :: igr_iter_solver - integer :: igr_order - logical :: igr - logical :: igr_pres_lim - logical :: viscous - #:endif + #:include 'generated_case_opt_decls.fpp' $:GPU_DECLARE(create='[int_comp, ic_eps, ic_beta]') integer :: hyper_model !< hyperelasticity solver algorithm @@ -256,32 +206,25 @@ module m_global_parameters ! END: Simulation Algorithm Parameters - ! Fluids Physical Parameters + ! Fluids Physical Parameters fluid_pp, bub_pp: auto-generated in generated_decls.fpp - type(physical_parameters), dimension(num_fluids_max) :: fluid_pp !< Stiffened gas EOS parameters and Reynolds numbers per fluid - ! Subgrid Bubble Parameters - type(subgrid_bubble_physical_parameters) :: bub_pp - integer :: fd_number !< Finite-difference half-stencil size: MAX(1, fd_order/2) + integer :: fd_number !< Finite-difference half-stencil size: MAX(1, fd_order/2) $:GPU_DECLARE(create='[fd_order, fd_number]') - type(vec3_dt), dimension(num_probes_max) :: probe - type(integral_parameters), dimension(num_probes_max) :: integral + ! probe, integral: auto-generated in generated_decls.fpp !> @name Reference density and pressure for Tait EOS !> @{ $:GPU_DECLARE(create='[rhoref, pref]') !> @name Immersed Boundaries + !> patch_ib, ib_airfoil, stl_models, particle_cloud: auto-generated in generated_decls.fpp !> @{ - type(ib_patch_parameters), dimension(num_ib_patches_max_namelist) :: patch_ib !< Immersed boundary patch parameters integer, dimension(num_local_ibs_max) :: local_ib_patch_ids !< lookup table of IBs in the local compute domain - type(particle_cloud_parameters), dimension(num_particle_clouds_max) :: particle_cloud !< Particle bed specifications integer, allocatable, dimension(:,:,:) :: ib_neighbor_ranks !< MPI ranks of neighborhood domains, indexed (-N:N,-N:N,-N:N) - type(ib_airfoil_parameters), dimension(num_ib_airfoils_max) :: ib_airfoil !< Per-airfoil NACA user inputs (namelist) type(ib_airfoil_grid), dimension(num_ib_airfoils_max) :: ib_airfoil_grids !< Per-airfoil computed surface grids - type(ib_stl_parameters), dimension(num_stl_models_max) :: stl_models !< Per-STL model parameters (namelist) - $:GPU_DECLARE(create='[ib, num_ibs, patch_ib, ib_airfoil, ib_airfoil_grids]') + $:GPU_DECLARE(create='[ib, num_ibs, ib_airfoil_grids]') $:GPU_DECLARE(create='[ib_coefficient_of_friction]') !> @} @@ -323,8 +266,7 @@ module m_global_parameters $:GPU_DECLARE(create='[mom_sp, mom_3d]') !> @} - type(chemistry_parameters) :: chem_params - $:GPU_DECLARE(create='[chem_params]') + ! chem_params: auto-generated in generated_decls.fpp !> @name Physical bubble parameters (see Ando 2010, Preston 2007) !> @{ @@ -343,11 +285,8 @@ module m_global_parameters $:GPU_DECLARE(create='[R0ref, p0ref, rho0ref, T0ref, ss, pv, vd, mu_l, mu_v, mu_g, gam_v, gam_g, M_v, M_g, cp_v, cp_g, R_v, R_g]') !> @} - !> @name Acoustic acoustic_source parameters - !> @{ - type(acoustic_parameters), dimension(num_probes_max) :: acoustic !< Acoustic source parameters - !> @} - $:GPU_DECLARE(create='[acoustic_source, acoustic, num_source]') + ! acoustic: auto-generated in generated_decls.fpp + $:GPU_DECLARE(create='[acoustic_source, num_source]') !> @name Surface tension parameters !> @{ @@ -365,9 +304,9 @@ module m_global_parameters $:GPU_DECLARE(create='[pb_ts, mv_ts]') !> @name lagrangian subgrid bubble parameters + !> lag_params: auto-generated in generated_decls.fpp !> @{! - type(bubbles_lagrange_parameters) :: lag_params !< Lagrange bubbles' parameters - $:GPU_DECLARE(create='[bubbles_lagrange, lag_params]') + $:GPU_DECLARE(create='[bubbles_lagrange]') !> @} $:GPU_DECLARE(create='[Bx0]') @@ -566,7 +505,7 @@ contains #:if not MFC_CASE_OPTIMIZATION nb = 1 - recon_type = WENO_TYPE + recon_type = recon_type_weno weno_order = dflt_int muscl_order = dflt_int muscl_lim = dflt_int @@ -777,14 +716,14 @@ contains #:if not MFC_CASE_OPTIMIZATION ! Determining the degree of the WENO polynomials - if (recon_type == WENO_TYPE) then + if (recon_type == recon_type_weno) then weno_polyn = (weno_order - 1)/2 if (teno) then weno_num_stencils = weno_order - 3 else weno_num_stencils = weno_polyn end if - else if (recon_type == MUSCL_TYPE) then + else if (recon_type == recon_type_muscl) then muscl_polyn = muscl_order end if $:GPU_UPDATE(device='[weno_polyn, muscl_polyn]') diff --git a/src/simulation/m_muscl.fpp b/src/simulation/m_muscl.fpp index d604648a76..11791cd309 100644 --- a/src/simulation/m_muscl.fpp +++ b/src/simulation/m_muscl.fpp @@ -10,6 +10,8 @@ module m_muscl use m_derived_types use m_global_parameters use m_variables_conversion + use m_constants, only: muscl_order_first_order, muscl_order_second_order, muscl_lim_unlimited, muscl_lim_minmod, & + & muscl_lim_mc, muscl_lim_van_albada, muscl_lim_van_leer, muscl_lim_superbee #ifdef MFC_OpenACC use openacc #endif @@ -99,7 +101,7 @@ contains $:GPU_UPDATE(device='[is1_muscl, is2_muscl, is3_muscl]') - if (muscl_order == 1) then + if (muscl_order == muscl_order_first_order) then if (muscl_dir == 1) then $:GPU_PARALLEL_LOOP(collapse=4) do i = 1, ubound(v_vf, 1) @@ -145,7 +147,7 @@ contains v_size = ubound(v_vf, 1) $:GPU_UPDATE(device='[v_size]') - if (muscl_order /= 1) then + if (muscl_order /= muscl_order_first_order) then $:GPU_PARALLEL_LOOP(private='[j, k, l, i]', collapse=4) do i = 1, v_size do l = idwbuff(3)%beg, idwbuff(3)%end @@ -159,7 +161,7 @@ contains $:END_GPU_PARALLEL_LOOP() end if - if (muscl_order == 2) then + if (muscl_order == muscl_order_second_order) then ! MUSCL Reconstruction #:for MUSCL_DIR, XYZ, STENCIL_VAR, COORDS, X_BND, Y_BND, Z_BND in & [(1, 'x', 'j', '{STENCIL_IDX}, k, l', 'is1_muscl', 'is2_muscl', 'is3_muscl'), & @@ -177,28 +179,28 @@ contains slopeR = v_rs_ws_muscl(${SF('')}$, i) - v_rs_ws_muscl(${SF(' - 1')}$, i) slope = 0._wp - if (muscl_lim == 0) then ! unlimited (central difference) + if (muscl_lim == muscl_lim_unlimited) then ! unlimited (central difference) slope = 5e-1_wp*(slopeL + slopeR) - else if (muscl_lim == 1) then ! minmod + else if (muscl_lim == muscl_lim_minmod) then ! minmod if (slopeL*slopeR > muscl_eps) then slope = min(abs(slopeL), abs(slopeR)) end if if (slopeL < 0._wp) slope = -slope - else if (muscl_lim == 2) then ! MC + else if (muscl_lim == muscl_lim_mc) then ! MC if (slopeL*slopeR > muscl_eps) then slope = min(2._wp*abs(slopeL), 2._wp*abs(slopeR)) slope = min(slope, 5e-1_wp*(abs(slopeL) + abs(slopeR))) end if if (slopeL < 0._wp) slope = -slope - else if (muscl_lim == 3) then ! Van Albada + else if (muscl_lim == muscl_lim_van_albada) then ! Van Albada if (slopeL*slopeR > muscl_eps) then slope = ((slopeL + slopeR)*slopeL*slopeR)/(slopeL**2._wp + slopeR**2._wp) end if - else if (muscl_lim == 4) then ! Van Leer + else if (muscl_lim == muscl_lim_van_leer) then ! Van Leer if (slopeL*slopeR > muscl_eps) then slope = 2._wp*slopeL*slopeR/(slopeL + slopeR) end if - else if (muscl_lim == 5) then ! SUPERBEE + else if (muscl_lim == muscl_lim_superbee) then ! SUPERBEE if (slopeL*slopeR > muscl_eps) then slope = -1._wp*min(-min(2._wp*abs(slopeL), abs(slopeR)), -min(abs(slopeL), & & 2._wp*abs(slopeR))) diff --git a/src/simulation/m_qbmm.fpp b/src/simulation/m_qbmm.fpp index 8465d4d349..abb77a1f35 100644 --- a/src/simulation/m_qbmm.fpp +++ b/src/simulation/m_qbmm.fpp @@ -14,6 +14,7 @@ module m_qbmm use m_variables_conversion use m_helper_basic use m_helper + use m_constants, only: bubble_model_keller_miksis, bubble_model_rayleigh_plesset implicit none @@ -43,10 +44,10 @@ contains integer :: i1, i2, q, i, j #:if not MFC_CASE_OPTIMIZATION - if (bubble_model == 2) then + if (bubble_model == bubble_model_keller_miksis) then ! Keller-Miksis without viscosity/surface tension nterms = 32 - else if (bubble_model == 3) then + else if (bubble_model == bubble_model_rayleigh_plesset) then ! Rayleigh-Plesset with viscosity/surface tension nterms = 7 end if @@ -64,7 +65,7 @@ contains do q = 1, nb do i1 = 0, 2; do i2 = 0, 2 if ((i1 + i2) <= 2) then - if (bubble_model == 3) then + if (bubble_model == bubble_model_rayleigh_plesset) then momrhs(1, i1, i2, 1, q) = -1._wp + i1 momrhs(2, i1, i2, 1, q) = -1._wp + i2 momrhs(3, i1, i2, 1, q) = 0._wp @@ -98,7 +99,7 @@ contains momrhs(1, i1, i2, 7, q) = -1._wp + i1 momrhs(2, i1, i2, 7, q) = -1._wp + i2 momrhs(3, i1, i2, 7, q) = 0._wp - else if (bubble_model == 2) then + else if (bubble_model == bubble_model_keller_miksis) then ! KM with approximation of 1/(1-V/C) = 1+V/C momrhs(1, i1, i2, 1, q) = -1._wp + i1 momrhs(2, i1, i2, 1, q) = 1._wp + i2 @@ -235,7 +236,7 @@ contains do q = 1, nb do i1 = 0, 2; do i2 = 0, 2 if ((i1 + i2) <= 2) then - if (bubble_model == 3) then + if (bubble_model == bubble_model_rayleigh_plesset) then momrhs(1, i1, i2, 1, q) = -1._wp + i1 momrhs(2, i1, i2, 1, q) = -1._wp + i2 momrhs(3, i1, i2, 1, q) = 0._wp @@ -269,7 +270,7 @@ contains momrhs(1, i1, i2, 7, q) = -1._wp + i1 momrhs(2, i1, i2, 7, q) = -1._wp + i2 momrhs(3, i1, i2, 7, q) = 0._wp - else if (bubble_model == 2) then + else if (bubble_model == bubble_model_keller_miksis) then ! KM with approximation of 1/(1-V/C) = 1+V/C momrhs(1, i1, i2, 1, q) = -1._wp + i1 momrhs(2, i1, i2, 1, q) = 1._wp + i2 @@ -599,7 +600,7 @@ contains do i2 = 0, 2; do i1 = 0, 2 if ((i1 + i2) <= 2) then - if (bubble_model == 3) then + if (bubble_model == bubble_model_rayleigh_plesset) then #:if not MFC_CASE_OPTIMIZATION or nterms > 1 ! RPE coeffs(1, i1, i2) = -1._wp*i2*pres/rho @@ -610,7 +611,7 @@ contains if (.not. f_is_default(Web)) coeffs(6, i1, i2) = -2._wp*i2/Web/rho coeffs(7, i1, i2) = 0._wp #:endif - else if (bubble_model == 2) then + else if (bubble_model == bubble_model_keller_miksis) then ! KM with approximation of 1/(1-V/C) = 1+V/C #:if not MFC_CASE_OPTIMIZATION or nterms > 7 coeffs(1, i1, i2) = -3._wp*i2/2._wp @@ -678,7 +679,7 @@ contains do i2 = 0, 2; do i1 = 0, 2 if ((i1 + i2) <= 2) then - if (bubble_model == 3) then + if (bubble_model == bubble_model_rayleigh_plesset) then ! RPE #:if not MFC_CASE_OPTIMIZATION or nterms > 1 coeffs(1, i1, i2) = -1._wp*i2*pres/rho @@ -689,7 +690,7 @@ contains if (.not. f_is_default(Web)) coeffs(6, i1, i2) = -2._wp*i2/Web/rho coeffs(7, i1, i2) = i2*pv/rho #:endif - else if (bubble_model == 2) then + else if (bubble_model == bubble_model_keller_miksis) then ! KM with approximation of 1/(1-V/C) = 1+V/C #:if not MFC_CASE_OPTIMIZATION or nterms > 7 coeffs(1, i1, i2) = -3._wp*i2/2._wp @@ -770,7 +771,7 @@ contains pres = q_prim_vf(eqn_idx%E)%sf(id1, id2, id3) rho = q_prim_vf(eqn_idx%cont%beg)%sf(id1, id2, id3) - if (bubble_model == 2) then + if (bubble_model == bubble_model_keller_miksis) then n_tait = 1._wp/gammas(1) + 1._wp B_tait = pi_infs(1)*(n_tait - 1)/n_tait c = n_tait*(pres + B_tait)*(1._wp - alf)/(rho) @@ -829,7 +830,7 @@ contains $:GPU_LOOP(parallelism='[seq]') do j = 1, nterms select case (bubble_model) - case (3) + case (bubble_model_rayleigh_plesset) if (j == 3) then momsum = momsum + coeff(j, i1, i2)*(R0(q)**momrhs(3, i1, i2, j, & & q))*f_quad2D(abscX(:,q), abscY(:,q), wght_pb(:,q), & @@ -839,7 +840,7 @@ contains & q))*f_quad2D(abscX(:,q), abscY(:,q), wght(:,q), & & momrhs(:,i1, i2, j, q)) end if - case (2) + case (bubble_model_keller_miksis) if ((j >= 7 .and. j <= 9) .or. (j >= 22 .and. j <= 23) .or. (j >= 10 & & .and. j <= 11) .or. (j == 26)) then momsum = momsum + coeff(j, i1, i2)*(R0(q)**momrhs(3, i1, i2, j, & diff --git a/src/simulation/m_rhs.fpp b/src/simulation/m_rhs.fpp index 6ba0f67dbd..1e30f135be 100644 --- a/src/simulation/m_rhs.fpp +++ b/src/simulation/m_rhs.fpp @@ -14,7 +14,8 @@ module m_rhs use m_mpi_proxy use m_variables_conversion use m_weno - use m_constants, only: riemann_solver_hll, riemann_solver_hlld, model_eqns_6eq + use m_constants, only: riemann_solver_hll, riemann_solver_hlld, model_eqns_6eq, int_comp_mthinc, recon_type_weno, & + & recon_type_muscl use m_muscl use m_riemann_solvers use m_cbc @@ -570,7 +571,7 @@ contains call nvtxEndRange end if - if (int_comp == 2 .and. n > 0) then + if (int_comp == int_comp_mthinc .and. n > 0) then call nvtxStartRange("RHS-COMPRESSION-NORMALS") call s_compute_mthinc_normals(q_prim_qp%vf) call nvtxEndRange @@ -1600,7 +1601,7 @@ contains integer :: recon_dir !< Coordinate direction of the reconstruction integer :: i, j, k, l - #:for SCHEME, TYPE in [('weno','WENO_TYPE'), ('muscl','MUSCL_TYPE')] + #:for SCHEME, TYPE in [('weno','recon_type_weno'), ('muscl','recon_type_muscl')] if (recon_type == ${TYPE}$) then ! Reconstruction in s1-direction if (norm_dir == 1) then @@ -1635,7 +1636,7 @@ contains integer :: i, j, k, l ! Reconstruction in s1-direction - #:for SCHEME, TYPE in [('weno','WENO_TYPE'), ('muscl', 'MUSCL_TYPE')] + #:for SCHEME, TYPE in [('weno','recon_type_weno'), ('muscl', 'recon_type_muscl')] if (recon_type == ${TYPE}$) then if (norm_dir == 1) then is1 = idwbuff(1); is2 = idwbuff(2); is3 = idwbuff(3) diff --git a/src/simulation/m_riemann_solvers.fpp b/src/simulation/m_riemann_solvers.fpp index fdb456bec4..f3dc56b7ec 100644 --- a/src/simulation/m_riemann_solvers.fpp +++ b/src/simulation/m_riemann_solvers.fpp @@ -16,7 +16,8 @@ module m_riemann_solvers use m_variables_conversion use m_bubbles use m_constants, only: riemann_solver_hll, riemann_solver_hllc, riemann_solver_hlld, riemann_solver_lax_friedrichs, & - & model_eqns_5eq, model_eqns_6eq, model_eqns_4eq + & model_eqns_5eq, model_eqns_6eq, model_eqns_4eq, avg_state_roe, avg_state_arithmetic, wave_speeds_direct, & + & wave_speeds_pressure use m_bubbles_EE use m_surface_tension use m_helper_basic @@ -483,7 +484,7 @@ contains end if ! Wave speed estimates (wave_speeds=1: direct, wave_speeds=2: pressure-based) - if (wave_speeds == 1) then + if (wave_speeds == wave_speeds_direct) then if (mhd) then ! MHD: use fast magnetosonic speed s_L = min(vel_L(dir_idx(1)) - c_fast%L, vel_R(dir_idx(1)) - c_fast%R) @@ -517,7 +518,7 @@ contains s_S = (pres_R - pres_L + rho_L*vel_L(dir_idx(1))*(s_L - vel_L(dir_idx(1))) & & - rho_R*vel_R(dir_idx(1))*(s_R - vel_R(dir_idx(1))))/(rho_L*(s_L - vel_L(dir_idx(1))) & & - rho_R*(s_R - vel_R(dir_idx(1)))) - else if (wave_speeds == 2) then + else if (wave_speeds == wave_speeds_pressure) then pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(dir_idx(1)) - vel_R(dir_idx(1)))) pres_SR = pres_SL @@ -1973,7 +1974,7 @@ contains end if ! COMPUTING THE DIRECT WAVE SPEEDS - if (wave_speeds == 1) then + if (wave_speeds == wave_speeds_direct) then if (elasticity) then ! Elastic wave speed, Rodriguez et al. JCP (2019) s_L = min(vel_L(dir_idx(1)) - sqrt(c_L*c_L + (((4._wp*G_L)/3._wp) + tau_e_L(dir_idx_tau(1) & @@ -1995,7 +1996,7 @@ contains & - rho_R*vel_R(dir_idx(1))*(s_R - vel_R(dir_idx(1))))/(rho_L*(s_L & & - vel_L(dir_idx(1))) - rho_R*(s_R - vel_R(dir_idx(1)))) end if - else if (wave_speeds == 2) then + else if (wave_speeds == wave_speeds_pressure) then pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(dir_idx(1)) - vel_R(dir_idx(1)))) pres_SR = pres_SL @@ -2277,14 +2278,14 @@ contains call s_compute_speed_of_sound(pres_R, rho_avg, gamma_avg, pi_inf_R, H_avg, alpha_R, vel_avg_rms, & & 0._wp, c_avg, qv_avg) - if (wave_speeds == 1) then + if (wave_speeds == wave_speeds_direct) then s_L = min(vel_L(dir_idx(1)) - c_L, vel_R(dir_idx(1)) - c_R) s_R = max(vel_R(dir_idx(1)) + c_R, vel_L(dir_idx(1)) + c_L) s_S = (pres_R - pres_L + rho_L*vel_L(dir_idx(1))*(s_L - vel_L(dir_idx(1))) & & - rho_R*vel_R(dir_idx(1))*(s_R - vel_R(dir_idx(1))))/(rho_L*(s_L - vel_L(dir_idx(1))) & & - rho_R*(s_R - vel_R(dir_idx(1)))) - else if (wave_speeds == 2) then + else if (wave_speeds == wave_speeds_pressure) then pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(dir_idx(1)) - vel_R(dir_idx(1)))) pres_SR = pres_SL @@ -2522,7 +2523,7 @@ contains H_L = (E_L + pres_L)/rho_L H_R = (E_R + pres_R)/rho_R - if (avg_state == 2) then + if (avg_state == avg_state_arithmetic) then $:GPU_LOOP(parallelism='[seq]') do i = 1, nb R0_L(i) = qL_prim_rsx_vf(${SF('')}$, rs(i)) @@ -2634,14 +2635,14 @@ contains @:compute_low_Mach_correction() end if - if (wave_speeds == 1) then + if (wave_speeds == wave_speeds_direct) then s_L = min(vel_L(dir_idx(1)) - c_L, vel_R(dir_idx(1)) - c_R) s_R = max(vel_R(dir_idx(1)) + c_R, vel_L(dir_idx(1)) + c_L) s_S = (pres_R - pres_L + rho_L*vel_L(dir_idx(1))*(s_L - vel_L(dir_idx(1))) & & - rho_R*vel_R(dir_idx(1))*(s_R - vel_R(dir_idx(1))))/(rho_L*(s_L - vel_L(dir_idx(1))) & & - rho_R*(s_R - vel_R(dir_idx(1)))) - else if (wave_speeds == 2) then + else if (wave_speeds == wave_speeds_pressure) then pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(dir_idx(1)) - vel_R(dir_idx(1)))) pres_SR = pres_SL @@ -2696,7 +2697,7 @@ contains ! Include p_tilde - if (avg_state == 2) then + if (avg_state == avg_state_arithmetic) then if (alpha_L(num_fluids) < small_alf .or. R3Lbar < small_alf) then pres_L = pres_L - alpha_L(num_fluids)*pres_L else @@ -3052,7 +3053,7 @@ contains @:compute_low_Mach_correction() end if - if (wave_speeds == 1) then + if (wave_speeds == wave_speeds_direct) then if (elasticity) then ! Elastic wave speed, Rodriguez et al. JCP (2019) s_L = min(vel_L(dir_idx(1)) - sqrt(c_L*c_L + (((4._wp*G_L)/3._wp) + tau_e_L(dir_idx_tau(1) & @@ -3074,7 +3075,7 @@ contains & - rho_R*vel_R(dir_idx(1))*(s_R - vel_R(dir_idx(1))))/(rho_L*(s_L & & - vel_L(dir_idx(1))) - rho_R*(s_R - vel_R(dir_idx(1)))) end if - else if (wave_speeds == 2) then + else if (wave_speeds == wave_speeds_pressure) then pres_SL = 5.e-1_wp*(pres_L + pres_R + rho_avg*c_avg*(vel_L(dir_idx(1)) - vel_R(dir_idx(1)))) pres_SR = pres_SL diff --git a/src/simulation/m_start_up.fpp b/src/simulation/m_start_up.fpp index 7a4a1a7df8..e0976f9270 100644 --- a/src/simulation/m_start_up.fpp +++ b/src/simulation/m_start_up.fpp @@ -51,7 +51,7 @@ module m_start_up use m_body_forces use m_sim_helpers use m_igr - use m_constants, only: model_eqns_6eq, time_stepper_rk1, time_stepper_rk2, time_stepper_rk3 + use m_constants, only: model_eqns_6eq, time_stepper_rk1, time_stepper_rk2, time_stepper_rk3, recon_type_weno, recon_type_muscl implicit none @@ -907,9 +907,9 @@ contains call s_initialize_igr_module() end if if (.not. igr) then - if (recon_type == WENO_TYPE) then + if (recon_type == recon_type_weno) then call s_initialize_weno_module() - else if (recon_type == MUSCL_TYPE) then + else if (recon_type == recon_type_muscl) then call s_initialize_muscl_module() end if call s_initialize_cbc_module() @@ -1081,9 +1081,9 @@ contains else call s_finalize_cbc_module() call s_finalize_riemann_solvers_module() - if (recon_type == WENO_TYPE) then + if (recon_type == recon_type_weno) then call s_finalize_weno_module() - else if (recon_type == MUSCL_TYPE) then + else if (recon_type == recon_type_muscl) then call s_finalize_muscl_module() end if end if diff --git a/src/simulation/m_thinc.fpp b/src/simulation/m_thinc.fpp index 3947df1032..5f9d1c7d54 100644 --- a/src/simulation/m_thinc.fpp +++ b/src/simulation/m_thinc.fpp @@ -16,6 +16,7 @@ module m_thinc use m_derived_types use m_global_parameters use m_helper + use m_constants, only: int_comp_mthinc #ifdef MFC_OpenACC use openacc @@ -194,7 +195,7 @@ contains subroutine s_initialize_thinc_module() - if (int_comp == 2) then + if (int_comp == int_comp_mthinc) then @:ALLOCATE(mthinc_nhat(1:3, idwbuff(1)%beg:idwbuff(1)%end, idwbuff(2)%beg:idwbuff(2)%end, & & idwbuff(3)%beg:idwbuff(3)%end)) @:ALLOCATE(mthinc_d(idwbuff(1)%beg:idwbuff(1)%end, idwbuff(2)%beg:idwbuff(2)%end, idwbuff(3)%beg:idwbuff(3)%end)) @@ -320,7 +321,7 @@ contains aCR = v_rs_ws(${SF(' + 1')}$, eqn_idx%adv%beg) if (aC >= ic_eps .and. aC <= 1._wp - ic_eps) then - if (int_comp == 2 .and. n > 0) then ! MTHINC + if (int_comp == int_comp_mthinc .and. n > 0) then ! MTHINC ! Map reshaped (j,k,l) to physical (ix,iy,iz) nh1 = mthinc_nhat(1, j, k, l) @@ -406,7 +407,7 @@ contains subroutine s_finalize_thinc_module() - if (int_comp == 2) then + if (int_comp == int_comp_mthinc) then @:DEALLOCATE(mthinc_nhat) @:DEALLOCATE(mthinc_d) end if diff --git a/src/simulation/m_viscous.fpp b/src/simulation/m_viscous.fpp index 4df940985f..d208428900 100644 --- a/src/simulation/m_viscous.fpp +++ b/src/simulation/m_viscous.fpp @@ -46,7 +46,7 @@ module m_viscous use m_muscl use m_helper use m_finite_differences - use m_constants, only: model_eqns_5eq + use m_constants, only: model_eqns_5eq, recon_type_weno, recon_type_muscl use m_hb_function private; public s_get_viscous, s_compute_viscous_stress_cylindrical_boundary, s_initialize_viscous_module, & @@ -848,7 +848,7 @@ contains integer :: recon_dir !< Coordinate direction of the WENO reconstruction integer :: i, j, k, l - #:for SCHEME, TYPE in [('weno','WENO_TYPE'), ('muscl','MUSCL_TYPE')] + #:for SCHEME, TYPE in [('weno','recon_type_weno'), ('muscl','recon_type_muscl')] if (recon_type == ${TYPE}$) then ! Reconstruction in s1-direction @@ -932,7 +932,7 @@ contains integer :: recon_dir !< Coordinate direction of the WENO reconstruction integer :: i, j, k, l - #:for SCHEME, TYPE in [('weno','WENO_TYPE'), ('muscl','MUSCL_TYPE')] + #:for SCHEME, TYPE in [('weno','recon_type_weno'), ('muscl','recon_type_muscl')] if (recon_type == ${TYPE}$) then ! Reconstruction in s1-direction diff --git a/toolchain/mfc/lint_docs.py b/toolchain/mfc/lint_docs.py index a5f88c7136..e852284821 100644 --- a/toolchain/mfc/lint_docs.py +++ b/toolchain/mfc/lint_docs.py @@ -67,30 +67,8 @@ # Hardcoded Fortran constants (not case-file params) "init_dir", "zeros_default", - # Enumerated value names used as prose examples (not parameter names) - "hll", - "hllc", - # Analytic expression language: Fortran intrinsics and module name (not case params) + # Analytic expression language: module name (not a case param) "m_constants", - "sin", - "cos", - "tan", - "asin", - "acos", - "atan", - "atan2", - "sinh", - "cosh", - "tanh", - "exp", - "log", - "log10", - "sqrt", - "abs", - "min", - "max", - "mod", - "sign", } # Docs to check for parameter references, with per-file skip sets @@ -196,11 +174,19 @@ def check_param_refs(repo_root: Path) -> list[str]: if toolchain_dir not in sys.path: sys.path.insert(0, toolchain_dir) try: + from mfc.analytic_expr import INTRINSICS, PASSTHROUGH from mfc.params import REGISTRY + from mfc.params.definitions import CONSTRAINTS except ImportError: print(" Warning: could not import REGISTRY, skipping parameter check") return [] + # Collect all enumerated value names from CONSTRAINTS (e.g. "hll", "hllc", "rk3") + _constraint_names: set[str] = set() + for _entry in CONSTRAINTS.values(): + if isinstance(_entry, dict) and "names" in _entry: + _constraint_names.update(_entry["names"].keys()) + valid_params = set(REGISTRY.all_params.keys()) # Build set of sub-parameter base names (strip trailing (N) indexes) _sub_raw = {p.split("%")[-1] for p in valid_params if "%" in p} @@ -234,6 +220,8 @@ def check_param_refs(repo_root: Path) -> list[str]: continue if param in extra_skip: continue + if param in _constraint_names or param in INTRINSICS or param in PASSTHROUGH: + continue if "(" in param or ")" in param: continue if "[" in param: diff --git a/toolchain/mfc/params/definitions.py b/toolchain/mfc/params/definitions.py index e9671ca5b8..5d0c0e55b0 100644 --- a/toolchain/mfc/params/definitions.py +++ b/toolchain/mfc/params/definitions.py @@ -1129,6 +1129,29 @@ def _init_registry(): "vel_wrt": "3", } +# Derived-type namelist variables whose Fortran declarations come from generated_decls.fpp. +# Maps variable name -> (fortran_type, dimension_expr_or_None, sim_gpu_declare, doxygen_desc_or_None). +# sim_gpu_declare=True means the generator emits this variable's $:GPU_DECLARE line +# alongside its declaration; the variable must then be REMOVED from any grouped +# GPU_DECLARE list in src/simulation/m_global_parameters.fpp (multiple declare +# directives accumulate identically in OpenACC and OpenMP). +TYPED_DECLS: dict[str, tuple] = { + "fluid_pp": ("type(physical_parameters)", "num_fluids_max", False, "Per-fluid stiffened-gas EOS parameters, Reynolds numbers, and shear modulus"), + "bub_pp": ("type(subgrid_bubble_physical_parameters)", None, False, "Subgrid bubble physical parameters"), + "patch_icpp": ("type(ic_patch_parameters)", "num_patches_max", False, "IC patch parameters"), + "patch_bc": ("type(bc_patch_parameters)", "num_bc_patches_max", False, "Boundary condition patch parameters"), + "patch_ib": ("type(ib_patch_parameters)", "num_ib_patches_max_namelist", True, "Immersed boundary patch parameters"), + "ib_airfoil": ("type(ib_airfoil_parameters)", "num_ib_airfoils_max", True, "Per-airfoil NACA user inputs"), + "stl_models": ("type(ib_stl_parameters)", "num_stl_models_max", False, "Per-STL model parameters"), + "probe": ("type(vec3_dt)", "num_probes_max", False, None), + "integral": ("type(integral_parameters)", "num_probes_max", False, None), + "acoustic": ("type(acoustic_parameters)", "num_probes_max", True, "Acoustic source parameters"), + "chem_params": ("type(chemistry_parameters)", None, True, None), + "lag_params": ("type(bubbles_lagrange_parameters)", None, True, "Lagrange bubbles' parameters"), + "particle_cloud": ("type(particle_cloud_parameters)", "num_particle_clouds_max", False, "Particle bed specifications"), + "simplex_params": ("type(simplex_noise_params)", None, False, None), +} + def _nv(targets: set, *names: str) -> None: for n in names: diff --git a/toolchain/mfc/params/generators/fortran_gen.py b/toolchain/mfc/params/generators/fortran_gen.py index c29a70d192..0ef8efe0c6 100644 --- a/toolchain/mfc/params/generators/fortran_gen.py +++ b/toolchain/mfc/params/generators/fortran_gen.py @@ -3,7 +3,7 @@ from pathlib import Path from typing import List, Tuple -from ..definitions import CASE_OPT_PARAMS, FORTRAN_ARRAY_DIMS, NAMELIST_VARS # noqa: F401 - triggers registry population +from ..definitions import CASE_OPT_PARAMS, FORTRAN_ARRAY_DIMS, NAMELIST_VARS, TYPED_DECLS # noqa: F401 - triggers registry population from ..registry import REGISTRY from ..schema import ParamDef, ParamType @@ -116,6 +116,9 @@ def generate_decls_fpp(target: str) -> str: continue if target == "sim" and name in CASE_OPT_PARAMS: continue + # TYPED_DECLS handles these; skip here to avoid double-emission. + if name in TYPED_DECLS: + continue if name in FORTRAN_ARRAY_DIMS: member = REGISTRY.all_params.get(f"{name}(1)") if member is None: @@ -130,6 +133,17 @@ def generate_decls_fpp(target: str) -> str: if any(k.startswith(f"{name}(") for k in REGISTRY.all_params): raise ValueError(f"{name!r} has indexed variants (e.g. {name}(1)) but is missing from " "FORTRAN_ARRAY_DIMS. Add it there with its Fortran dimension expression.") lines.append(f"{fortran_type_decl(param).ljust(_DECL_COL)}:: {name}") + for name, (ftype, dim, gpu, desc) in TYPED_DECLS.items(): + if name not in NAMELIST_VARS or target not in NAMELIST_VARS[name]: + continue + decl = f"{ftype}, dimension({dim})" if dim else ftype + padded = decl.ljust(_ARRAY_DECL_COL) + if not padded.endswith(" "): + padded += " " + doc = f" !< {desc}" if desc else "" + lines.append(f"{padded}:: {name}{doc}") + if gpu and target == "sim": + lines.append(f"$:GPU_DECLARE(create='[{name}]')") return "\n".join(lines) + "\n" @@ -147,6 +161,74 @@ def generate_constants_fpp() -> str: return "\n".join(lines) + "\n" +# case.py-computed extras that are not CASE_OPT_PARAMS but appear in the +# case-optimization declaration block (injected as Fypp #:set variables +# by toolchain/mfc/case.py, not from the namelist registry). +# Order matches the reference manual block. +CASE_OPT_EXTRA_LINES = [ + ("num_dims", "integer", "Number of spatial dimensions"), + ("num_vels", "integer", "Number of velocity components (different from num_dims for mhd)"), + ("weno_polyn", "integer", "Degree of the WENO polynomials"), + ("muscl_polyn", "integer", "Degree of the MUSCL polynomials"), + ("weno_num_stencils", "integer", "Number of stencils for WENO reconstruction"), + ("wenojs", "logical", "WENO-JS (default)"), +] + +_CASE_OPT_DECL_COL = 24 # '::' alignment for case-opt declarations + + +def generate_case_opt_decls_fpp() -> str: + """Return the case-optimization declaration block for src/simulation/m_global_parameters.fpp. + + Emits one #:if MFC_CASE_OPTIMIZATION / #:else / #:endif block containing: + - The CASE_OPT_EXTRA_LINES (case.py-computed: num_dims, num_vels, weno_polyn, + muscl_polyn, weno_num_stencils, wenojs) — stored as a literal list because + they are not namelist-registry parameters. + - Every CASE_OPT_PARAMS entry (excluding nb, which lives in a separate block) + driven by the registry for its Fortran type. + """ + params_to_emit = sorted(CASE_OPT_PARAMS - {"nb"}) + + if_lines = [] + else_lines = [] + + # Emit extras first (case.py-computed, not in registry) + for name, ftype, desc in CASE_OPT_EXTRA_LINES: + if ftype == "logical": + rhs = f"(${{{name}}}$ /= 0)" + if_lines.append(f" {ftype}, parameter :: {name} = {rhs} !< {desc}") + else: + if_lines.append(f" {ftype}, parameter :: {name} = ${{{name}}}$ !< {desc}") + else_lines.append(f" {ftype.ljust(_CASE_OPT_DECL_COL)}:: {name}") + + # Emit registry-driven CASE_OPT_PARAMS + for name in params_to_emit: + pdef = REGISTRY.get_param_def(name) + if pdef is None: + raise ValueError(f"CASE_OPT_PARAMS entry {name!r} not found in registry") + desc = pdef.description or "" + if pdef.param_type in _REAL_TYPES: + ftype = "real(wp)" + if_lines.append(f" {ftype}, parameter :: {name} = ${{{name}}}$ !< {desc}") + else_lines.append(f" {ftype.ljust(_CASE_OPT_DECL_COL)}:: {name}") + elif pdef.param_type == ParamType.LOG: + ftype = "logical" + rhs = f"(${{{name}}}$ /= 0)" + if_lines.append(f" {ftype}, parameter :: {name} = {rhs} !< {desc}") + else_lines.append(f" {ftype.ljust(_CASE_OPT_DECL_COL)}:: {name}") + else: + ftype = "integer" + if_lines.append(f" {ftype}, parameter :: {name} = ${{{name}}}$ !< {desc}") + else_lines.append(f" {ftype.ljust(_CASE_OPT_DECL_COL)}:: {name}") + + parts = [_HEADER.rstrip(), "#:if MFC_CASE_OPTIMIZATION"] + parts.extend(if_lines) + parts.append("#:else") + parts.extend(else_lines) + parts.append("#:endif") + return "\n".join(parts) + "\n" + + def resolve_namelist_content(fpp_path: Path) -> str: """Return the namelist content for an fpp file. @@ -163,10 +245,11 @@ def resolve_namelist_content(fpp_path: Path) -> str: def get_generated_files(build_dir: Path) -> List[Tuple[Path, str]]: - """Return (path, content) for all 9 generated .fpp files under build_dir. + """Return (path, content) for all 10 generated .fpp files under build_dir. Paths match the cmake include directory structure: build_dir/include/{full_target}/generated_{namelist,decls,constants}.fpp + The simulation target also gets generated_case_opt_decls.fpp. """ result = [] for short, full in TARGETS: @@ -174,4 +257,6 @@ def get_generated_files(build_dir: Path) -> List[Tuple[Path, str]]: result.append((inc / "generated_namelist.fpp", generate_namelist_fpp(short))) result.append((inc / "generated_decls.fpp", generate_decls_fpp(short))) result.append((inc / "generated_constants.fpp", generate_constants_fpp())) + sim_inc = build_dir / "include" / "simulation" + result.append((sim_inc / "generated_case_opt_decls.fpp", generate_case_opt_decls_fpp())) return result diff --git a/toolchain/mfc/params_tests/test_fortran_gen.py b/toolchain/mfc/params_tests/test_fortran_gen.py index 6546f8b425..6f007c0ead 100644 --- a/toolchain/mfc/params_tests/test_fortran_gen.py +++ b/toolchain/mfc/params_tests/test_fortran_gen.py @@ -140,14 +140,65 @@ def test_decls_array_dims(): assert "dimension(num_fluids_max)" in post and ":: alpha_wrt" in post assert "dimension(num_fluids_max)" in post and ":: alpha_rho_wrt" in post assert "dimension(3)" in post and ":: mom_wrt" in post - # Structs and families must NOT appear as bare scalar declarations - assert ":: fluid_pp" not in post + # fluid_pp is now emitted via TYPED_DECLS; bc_x still must not appear + assert ":: fluid_pp" in post assert ":: bc_x" not in post pre = generate_decls_fpp("pre") assert "dimension(num_fluids_max)" in pre and ":: fluid_rho" in pre +def test_generate_decls_emits_typed_declarations(): + from mfc.params.generators.fortran_gen import generate_decls_fpp + + pre = generate_decls_fpp("pre") + sim = generate_decls_fpp("sim") + post = generate_decls_fpp("post") + + # patch_icpp is pre-only + assert "patch_icpp" in pre and "patch_icpp" not in post and "patch_icpp" not in sim + + # fluid_pp appears in all three targets + for out in (pre, sim, post): + assert "fluid_pp" in out + + # acoustic is sim-only + assert "acoustic" in sim and "acoustic" not in pre + + # Exact line shape (harvested from the manual declaration it replaces). + # Note: _ARRAY_DECL_COL=36; the type+dim string is longer so no padding space before :: + assert "type(physical_parameters), dimension(num_fluids_max) :: fluid_pp" in sim + + # chem_params is sim-only with GPU declare + assert "chem_params" in sim and "chem_params" not in pre and "chem_params" not in post + assert "$:GPU_DECLARE(create='[chem_params]')" in sim + + # bub_pp is in all three targets (derived-type scalar, dim=None) + for out in (pre, sim, post): + assert "bub_pp" in out + + # lag_params is sim-only; the generator owns its GPU_DECLARE (removed from the + # grouped list in m_global_parameters.fpp by the Task-2 Fortran edit) + assert "lag_params" in sim and "lag_params" not in pre + for gpu_var in ("patch_ib", "ib_airfoil", "acoustic", "lag_params", "chem_params"): + assert f"$:GPU_DECLARE(create='[{gpu_var}]')" in sim + + # Doxygen descriptions from TYPED_DECLS are emitted as inline comments + assert ":: fluid_pp !< Per-fluid stiffened-gas EOS parameters, Reynolds numbers, and shear modulus" in sim + + # simplex_params is pre-only + assert "simplex_params" in pre and "simplex_params" not in sim + + # No TYPED_DECLS name emitted twice + for name in ("fluid_pp", "bub_pp", "chem_params", "lag_params"): + assert sim.count(f":: {name}") == 1, f"{name!r} emitted more than once in sim" + + # TYPED_DECLS keys must not overlap with FORTRAN_ARRAY_DIMS + from mfc.params.definitions import FORTRAN_ARRAY_DIMS, TYPED_DECLS + + assert not (set(TYPED_DECLS) & set(FORTRAN_ARRAY_DIMS)), "TYPED_DECLS and FORTRAN_ARRAY_DIMS share keys" + + def test_check_target_raises_on_bad_target(): import pytest @@ -159,17 +210,18 @@ def test_check_target_raises_on_bad_target(): generate_decls_fpp("bad") -def test_get_generated_files_returns_nine(): +def test_get_generated_files_returns_ten(): from pathlib import Path from mfc.params.generators.fortran_gen import get_generated_files files = get_generated_files(Path("/build")) - assert len(files) == 9 + assert len(files) == 10 paths = [str(p) for p, _ in files] assert any("pre_process/generated_namelist.fpp" in p for p in paths) assert any("simulation/generated_decls.fpp" in p for p in paths) assert any("post_process/generated_namelist.fpp" in p for p in paths) + assert any("simulation/generated_case_opt_decls.fpp" in p for p in paths) def test_generate_constants_fpp_content(): @@ -191,5 +243,47 @@ def test_get_generated_files_includes_constants(): files = get_generated_files(Path("/tmp/x")) names = {p.name for p, _ in files} - assert names == {"generated_namelist.fpp", "generated_decls.fpp", "generated_constants.fpp"} - assert len(files) == 9 + assert names == {"generated_namelist.fpp", "generated_decls.fpp", "generated_constants.fpp", "generated_case_opt_decls.fpp"} + assert len(files) == 10 + + +def test_generate_case_opt_decls_fpp(): + from mfc.params.generators.fortran_gen import generate_case_opt_decls_fpp + + out = generate_case_opt_decls_fpp() + + # Must have the case-opt guard structure + assert "#:if MFC_CASE_OPTIMIZATION" in out + assert "#:else" in out + assert "#:endif" in out + + # Representative registry-driven parameter lines (#:if branch) + assert "integer, parameter :: weno_order = ${weno_order}$" in out + assert "integer, parameter :: num_fluids = ${num_fluids}$" in out + assert "logical, parameter :: mapped_weno = (${mapped_weno}$ /= 0)" in out + assert "real(wp), parameter :: wenoz_q = ${wenoz_q}$" in out + + # #:else branch variable forms + assert "integer :: weno_order" in out + assert "logical :: mapped_weno" in out + assert "real(wp) :: wenoz_q" in out + + # case.py-computed extras + assert "integer, parameter :: num_dims = ${num_dims}$" in out + assert "integer, parameter :: weno_polyn = ${weno_polyn}$" in out + assert "logical, parameter :: wenojs = (${wenojs}$ /= 0)" in out + assert "integer :: num_dims" in out + assert "integer :: weno_polyn" in out + assert "logical :: wenojs" in out + + # nb must NOT be in this block (it lives in a separate block) + assert ":: nb" not in out + + # All CASE_OPT_PARAMS (minus nb) must appear in both branches + from mfc.params.definitions import CASE_OPT_PARAMS + + for name in CASE_OPT_PARAMS - {"nb"}: + assert f":: {name}" in out, f"{name!r} missing from generated_case_opt_decls" + + # AUTO-GENERATED header present + assert "AUTO-GENERATED" in out