diff --git a/applications/skeleton/source/algorithm/skeleton_constants_mod.x90 b/applications/skeleton/source/algorithm/skeleton_constants_mod.x90 index 4717f9b7c..0a7a0d82f 100644 --- a/applications/skeleton/source/algorithm/skeleton_constants_mod.x90 +++ b/applications/skeleton/source/algorithm/skeleton_constants_mod.x90 @@ -15,7 +15,7 @@ module skeleton_constants_mod ! Infrastructure - use constants_mod, only: str_def, i_def, l_def + use constants_mod, only: str_def, r_def, i_def, l_def use driver_modeldb_mod, only: modeldb_type use field_mod, only: field_type use fs_continuity_mod, only: W0, W1, W2, W2broken, & @@ -31,6 +31,7 @@ module skeleton_constants_mod use quadrature_rule_gaussian_mod, only: quadrature_rule_gaussian_type use timing_mod, only: start_timing, stop_timing, & tik, LPROF + use sci_mesh_enums_mod, only: get_mesh_enums ! Kernels use sci_compute_derham_matrices_kernel_mod, only: compute_derham_matrices_kernel_type @@ -109,6 +110,10 @@ contains integer(i_def) :: order_h, order_v integer(i_def) :: nqp_h_exact, nqp_v_exact + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius integer(tik) :: id @@ -116,10 +121,14 @@ contains call log_event( "Skeleton: creating runtime constants", LOG_LEVEL_INFO ) - order_h = modeldb%config%finite_element%element_order_h() - order_v = modeldb%config%finite_element%element_order_v() - nqp_h_exact = modeldb%config%finite_element%nqp_h_exact() - nqp_v_exact = modeldb%config%finite_element%nqp_v_exact() + call get_mesh_enums(mesh, geometry, topology) + + order_h = modeldb%config%finite_element%element_order_h() + order_v = modeldb%config%finite_element%element_order_v() + nqp_h_exact = modeldb%config%finite_element%nqp_h_exact() + nqp_v_exact = modeldb%config%finite_element%nqp_v_exact() + coord_system = modeldb%config%finite_element%coord_system() + scaled_radius = modeldb%config%planet%scaled_radius() !=========== Create function spaces required for setup ==================! @@ -167,13 +176,16 @@ contains call curl_inventory%add_operator( curl, w2_fs, w1_fs, mesh ) call broken_div_inventory%add_operator( broken_div, w3_fs, w2b_fs, mesh ) - call invoke( name = "create_de_rham_matrices", & - compute_derham_matrices_kernel_type(mm_w0, mm_w1, & - mm_w2, mm_w2b, & - mm_w3, mm_wtheta, & - grad, curl, & - div, broken_div, & - chi, panel_id, qr) ) + call invoke( name = "create_de_rham_matrices", & + compute_derham_matrices_kernel_type(mm_w0, mm_w1, & + mm_w2, mm_w2b, & + mm_w3, mm_wtheta, & + grad, curl, & + div, broken_div, & + chi, panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, qr) ) !===================== Add operators to inventories =======================! diff --git a/components/driver/source/driver_coordinates_mod.F90 b/components/driver/source/driver_coordinates_mod.F90 index d5c5d98a4..783acd40a 100644 --- a/components/driver/source/driver_coordinates_mod.F90 +++ b/components/driver/source/driver_coordinates_mod.F90 @@ -17,12 +17,13 @@ module driver_coordinates_mod xyz2alphabetar, alphabetar2xyz, & schmidt_transform_xyz, & inverse_schmidt_transform_xyz - - ! Configuration modules - use base_mesh_config_mod, only: geometry_planar, & + use sci_mesh_enums_mod, only: get_mesh_enums, & + geometry_planar, & geometry_spherical, & topology_fully_periodic, & topology_non_periodic + + ! Configuration modules use finite_element_config_mod, only: coord_system_xyz implicit none @@ -108,17 +109,7 @@ subroutine assign_coordinate_field(config, mesh, chi, panel_id) nullify( map, map_pid, dof_coords, reference_element ) - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() diff --git a/components/driver/source/driver_fem_mod.f90 b/components/driver/source/driver_fem_mod.f90 index 8d2bfe5e4..19441f458 100644 --- a/components/driver/source/driver_fem_mod.f90 +++ b/components/driver/source/driver_fem_mod.f90 @@ -25,6 +25,11 @@ module driver_fem_mod use mesh_collection_mod, only: mesh_collection use sci_chi_transform_mod, only: init_chi_transforms, & final_chi_transforms + use sci_mesh_enums_mod, only: get_mesh_enums, & + geometry_spherical, & + geometry_planar, & + topology_non_periodic, & + topology_fully_periodic ! Object types use config_mod, only: config_type @@ -33,10 +38,6 @@ module driver_fem_mod use inventory_by_mesh_mod, only: inventory_by_mesh_type ! Configuration modules - use base_mesh_config_mod, only: geometry_spherical, & - geometry_planar, & - topology_non_periodic, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_xyz, & coord_space_W0, & coord_space_Wchi, & @@ -121,17 +122,7 @@ subroutine init_fem(config, chi_inventory, panel_id_inventory) mesh => mesh_collection%get_mesh(all_mesh_names(i)) mesh_name = mesh%get_mesh_name() - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) ! Initialise coordinate transformations call init_chi_transforms( geometry, topology, & diff --git a/components/driver/source/driver_mesh_mod.f90 b/components/driver/source/driver_mesh_mod.f90 index 40157459a..fecb0497d 100644 --- a/components/driver/source/driver_mesh_mod.f90 +++ b/components/driver/source/driver_mesh_mod.f90 @@ -51,11 +51,11 @@ module driver_mesh_mod use global_mesh_collection_mod, only: global_mesh_collection use local_mesh_collection_mod, only: local_mesh_collection + use sci_mesh_enums_mod, only: geometry_spherical, & + topology_fully_periodic ! Configuration modules use finite_element_config_mod, only: cellshape_quadrilateral - use base_mesh_config_mod, only: geometry_spherical, & - topology_fully_periodic implicit none diff --git a/components/driver/unit-test/assign_coordinate_xyz_mod_test.pf b/components/driver/unit-test/assign_coordinate_xyz_mod_test.pf index 37eaac9ef..1ed944e46 100644 --- a/components/driver/unit-test/assign_coordinate_xyz_mod_test.pf +++ b/components/driver/unit-test/assign_coordinate_xyz_mod_test.pf @@ -9,7 +9,7 @@ module assign_coordinate_xyz_mod_test use constants_mod, only : r_def, i_def - use base_mesh_config_mod, only: geometry_planar, topology_fully_periodic + use sci_mesh_enums_mod, only: geometry_planar, topology_fully_periodic use funit implicit none diff --git a/components/science/Makefile b/components/science/Makefile index ac6023dd5..52826618a 100644 --- a/components/science/Makefile +++ b/components/science/Makefile @@ -71,9 +71,7 @@ document-api: api-documentation # unit-tests/%: export BIN_DIR ?= $(PROJECT_DIR)/test unit-tests/%: export CXX_LINK = TRUE -#unit-tests/%: export EXTERNAL_STATIC_LIBRARIES += pfunit unit-tests/%: export IMPORT_PARTS = $(CORE_ROOT_DIR)/infrastructure \ -# $(CORE_ROOT_DIR)/components/lfric-xios \ $(CORE_ROOT_DIR)/components/inventory \ $(CORE_ROOT_DIR)/components/driver unit-tests/%: export META_FILE_DIR = rose-meta/lfric-science/HEAD diff --git a/components/science/source/algorithm/sci_fem_constants_mod.x90 b/components/science/source/algorithm/sci_fem_constants_mod.x90 index 28c31d995..bafc9ffe8 100644 --- a/components/science/source/algorithm/sci_fem_constants_mod.x90 +++ b/components/science/source/algorithm/sci_fem_constants_mod.x90 @@ -15,6 +15,7 @@ module sci_fem_constants_mod ! Infrastructure + use config_mod, only: config_type use constants_mod, only: i_def, r_def, r_second, & str_def, l_def use field_mod, only: field_type @@ -31,15 +32,11 @@ module sci_fem_constants_mod use quadrature_rule_gaussian_mod, only: quadrature_rule_gaussian_type use timing_mod, only: start_timing, stop_timing, & tik, LPROF - - ! Configuration - use finite_element_config_mod, only: element_order_h, & - element_order_v, & - nqp_h_exact, & - nqp_v_exact + use sci_mesh_enums_mod, only: get_mesh_enums, & + topology_fully_periodic ! Other algorithms - use sci_geometric_constants_mod, only: get_coordinates, & + use sci_geometric_constants_mod, only: get_coordinates, & get_panel_id ! Kernels use sci_compute_div_operator_kernel_mod, & @@ -48,6 +45,12 @@ module sci_fem_constants_mod use dg_inc_matrix_vector_kernel_mod, only: dg_inc_matrix_vector_kernel_type use sci_multiplicity_kernel_mod, only: multiplicity_kernel_type + ! Configuration modules + use finite_element_config_mod, only: element_order_h, & + element_order_v, & + nqp_h_exact, & + nqp_v_exact + implicit none private @@ -197,10 +200,11 @@ contains end function get_qr_fv !> @brief Returns a pointer to a finite element mass matrix + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @param[in] space The space of the desired mass matrix - !> @param[in] mesh_id The ID of the mesh to get the object for !> @return The mass matrix operator - function get_mass_matrix_fe(space, mesh_id) result(mm_op) + function get_mass_matrix_fe(config, mesh, space) result(mm_op) use sci_compute_mass_matrix_kernel_w_scalar_mod, & only: compute_mass_matrix_kernel_w_scalar_type @@ -210,13 +214,13 @@ contains only: compute_mass_matrix_kernel_w1_type use sci_edge_lump_w2_mass_matrix_kernel_mod, & only: edge_lump_w2_mass_matrix_kernel_type - use base_mesh_config_mod, only: topology, topology_fully_periodic implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(operator_type), pointer :: mm_op logical(kind=l_def) :: constant_exists type(field_type), pointer :: chi(:) @@ -230,9 +234,19 @@ contains integer(kind=i_def), parameter :: stencil_depth = 1_i_def integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + ! If running at lowest order, use finite volume if (element_order_h == 0 .and. element_order_v == 0) then - mm_op => get_mass_matrix_fv(space, mesh_id) + mm_op => get_mass_matrix_fv(config, mesh, space) return end if @@ -266,8 +280,7 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist @@ -285,16 +298,28 @@ contains case (W1) call invoke( name="compute_w1_mass_matrix_fe", & compute_mass_matrix_kernel_w1_type(mm_op, chi, & - panel_id, qr_ptr) ) + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + qr_ptr) ) case (W2, W2broken, W2H, W2V) call invoke( name="compute_w2_mass_matrix_fe", & compute_mass_matrix_kernel_w2_type(mm_op, chi, & - panel_id, qr_ptr) ) + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + qr_ptr) ) case (W3, Wtheta) call invoke( name= "compute_scalar_mass_matrix_fe", & compute_mass_matrix_kernel_w_scalar_type(mm_op, chi, & panel_id, & extend_mesh, & + geometry, & + topology, & + coord_system, & + scaled_radius, & qr_ptr) ) end select ! Lump W2 mass matrix along the edge of the domain to avoid boundary @@ -315,10 +340,11 @@ contains end function get_mass_matrix_fe !> @brief Returns a pointer to a finite volume mass matrix + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @param[in] space The space of the desired mass matrix - !> @param[in] mesh_id The ID of the mesh to get the object for !> @return The mass matrix operator - function get_mass_matrix_fv(space, mesh_id) result(mm_op) + function get_mass_matrix_fv(config, mesh, space) result(mm_op) use sci_compute_mass_matrix_kernel_w_scalar_mod, & only: compute_mass_matrix_kernel_w_scalar_type @@ -328,13 +354,13 @@ contains only: compute_mass_matrix_kernel_w1_type use sci_edge_lump_w2_mass_matrix_kernel_mod, & only: edge_lump_w2_mass_matrix_kernel_type - use base_mesh_config_mod, only: topology, topology_fully_periodic implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(operator_type), pointer :: mm_op logical(kind=l_def) :: constant_exists type(field_type), pointer :: chi(:) @@ -348,6 +374,16 @@ contains integer(kind=i_def), parameter :: stencil_depth = 1_i_def integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + ! Point to appropriate inventory for this space select case (space) case (W1) @@ -378,8 +414,7 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist @@ -396,16 +431,28 @@ contains case (W1) call invoke( name="compute_w1_mass_matrix_fv", & compute_mass_matrix_kernel_w1_type(mm_op, chi, & - panel_id, qr_ptr) ) + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + qr_ptr) ) case (W2, W2broken, W2H, W2V) call invoke( name="compute_w2_mass_matrix_fv", & compute_mass_matrix_kernel_w2_type(mm_op, chi, & - panel_id, qr_ptr) ) + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + qr_ptr) ) case (W3, Wtheta) call invoke( name= "compute_scalar_mass_matrix_fv", & compute_mass_matrix_kernel_w_scalar_type(mm_op, chi, & panel_id, & extend_mesh, & + geometry, & + topology, & + coord_system, & + scaled_radius, & qr_ptr) ) end select ! Lump W2 mass matrix along the edge of the domain to avoid boundary @@ -426,18 +473,21 @@ contains end function get_mass_matrix_fv !> @brief Returns a pointer to a diagonal finite element mass matrix + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @param[in] space The space of the desired diagonal mass matrix - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @return The mass matrix diagonal field - function get_mass_matrix_diagonal_fe(space, mesh_id) result(diagonal_mm) + function get_mass_matrix_diagonal_fe(config, mesh, space) result(diagonal_mm) - use sci_mm_diagonal_kernel_mod, only: mm_diagonal_kernel_type + use sci_mm_diagonal_kernel_mod, only: mm_diagonal_kernel_type implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(field_type), pointer :: diagonal_mm type(operator_type), pointer :: mass_matrix logical(kind=l_def) :: constant_exists @@ -448,7 +498,7 @@ contains ! If running at lowest order, use finite volume if (element_order_h == 0 .and. element_order_v == 0) then - diagonal_mm => get_mass_matrix_diagonal_fv(space, mesh_id) + diagonal_mm => get_mass_matrix_diagonal_fv(config, mesh, space) return end if @@ -482,19 +532,18 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist - mass_matrix => get_mass_matrix_fe(space, mesh_id) + mass_matrix => get_mass_matrix_fe(config, mesh, space) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) fs => function_space_collection%get_fs(mesh, element_order_h, & element_order_v, space) call inventory%add_field(diagonal_mm, fs, mesh) - call invoke ( setval_c(diagonal_mm, 0.0_r_def), & + call invoke ( setval_c(diagonal_mm, 0.0_r_def), & mm_diagonal_kernel_type(diagonal_mm, mass_matrix) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) @@ -506,18 +555,20 @@ contains end function get_mass_matrix_diagonal_fe !> @brief Returns a pointer to a diagonal finite volume mass matrix + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @param[in] space The space of the desired diagonal mass matrix - !> @param[in] mesh_id The ID of the mesh to get the object for !> @return The mass matrix diagonal field - function get_mass_matrix_diagonal_fv(space, mesh_id) result(diagonal_mm) + function get_mass_matrix_diagonal_fv(config, mesh, space) result(diagonal_mm) use sci_mm_diagonal_kernel_mod, only: mm_diagonal_kernel_type implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(field_type), pointer :: diagonal_mm type(operator_type), pointer :: mass_matrix logical(kind=l_def) :: constant_exists @@ -556,18 +607,17 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist - mass_matrix => get_mass_matrix_fv(space, mesh_id) + mass_matrix => get_mass_matrix_fv(config, mesh, space) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) fs => function_space_collection%get_fs(mesh, 0, 0, space) call inventory%add_field(diagonal_mm, fs, mesh) - call invoke ( setval_c(diagonal_mm, 0.0_r_def), & + call invoke ( setval_c(diagonal_mm, 0.0_r_def), & mm_diagonal_kernel_type(diagonal_mm, mass_matrix) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) @@ -580,18 +630,20 @@ contains !> @brief Returns a pointer to a field for a finite element inverse lumped !! mass matrix + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @param[in] space The space of the desired lumped inverse mass matrix - !> @param[in] mesh_id The ID of the mesh to get the object for !> @return The lumped inverse mass matrix - function get_inverse_lumped_mass_matrix_fe(space, mesh_id) result(mm_linv) + function get_inverse_lumped_mass_matrix_fe(config, mesh, space) result(mm_linv) use dg_matrix_vector_kernel_mod, only: dg_matrix_vector_kernel_type implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(field_type) :: ones type(field_type), pointer :: mm_linv type(operator_type), pointer :: mass_matrix @@ -604,7 +656,7 @@ contains ! If running at lowest order, use finite volume if (element_order_h == 0 .and. element_order_v == 0) then - mm_linv => get_inverse_lumped_mass_matrix_fv(space, mesh_id) + mm_linv => get_inverse_lumped_mass_matrix_fv(config, mesh, space) return end if @@ -626,12 +678,11 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist - mass_matrix => get_mass_matrix_fe(space, mesh_id) + mass_matrix => get_mass_matrix_fe(config, mesh, space) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) @@ -664,18 +715,20 @@ contains !> @brief Returns a pointer to a field for a finite volume inverse lumped !! mass matrix + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @param[in] space The space of the desired lumped inverse mass matrix - !> @param[in] mesh_id The ID of the mesh to get the object for !> @return The lumped inverse mass matrix - function get_inverse_lumped_mass_matrix_fv(space, mesh_id) result(mm_linv) + function get_inverse_lumped_mass_matrix_fv(config, mesh, space) result(mm_linv) use dg_matrix_vector_kernel_mod, only: dg_matrix_vector_kernel_type implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(field_type) :: ones type(field_type), pointer :: mm_linv type(operator_type), pointer :: mass_matrix @@ -704,12 +757,11 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist - mass_matrix => get_mass_matrix_fv(space, mesh_id) + mass_matrix => get_mass_matrix_fv(config, mesh, space) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) @@ -741,17 +793,20 @@ contains !> @brief Returns a pointer to the finite element inverse W3 mass matrix !! operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for + !> @param[in] space The space of the desired inverse mass matrix !> @return The inverse W3 mass matrix operator - function get_inverse_mass_matrix_fe(space, mesh_id) result(mm_inv) + function get_inverse_mass_matrix_fe(config, mesh, space) result(mm_inv) use sci_invert_local_operator_kernel_mod, only: invert_local_operator_kernel_type implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(operator_type), pointer :: mass_matrix type(operator_type), pointer :: mm_inv logical(kind=l_def) :: constant_exists @@ -762,7 +817,7 @@ contains ! If running at lowest order, use finite volume if (element_order_h == 0 .and. element_order_v == 0) then - mm_inv => get_inverse_mass_matrix_fv(space, mesh_id) + mm_inv => get_inverse_mass_matrix_fv(config, mesh, space) return end if @@ -784,12 +839,11 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist - mass_matrix => get_mass_matrix_fe(space, mesh_id) + mass_matrix => get_mass_matrix_fe(config, mesh, space) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) @@ -810,17 +864,20 @@ contains !> @brief Returns a pointer to the finite volume inverse W3 mass matrix !! operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for + !> @param[in] space The space of the desired operator !> @return The inverse W3 mass matrix operator - function get_inverse_mass_matrix_fv(space, mesh_id) result(mm_inv) + function get_inverse_mass_matrix_fv(config, mesh, space) result(mm_inv) use sci_invert_local_operator_kernel_mod, only: invert_local_operator_kernel_type implicit none - integer(kind=i_def), intent(in) :: space - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + integer(kind=i_def), intent(in) :: space + type(operator_type), pointer :: mass_matrix type(operator_type), pointer :: mm_inv logical(kind=l_def) :: constant_exists @@ -847,12 +904,11 @@ contains call inventory%initialise(name=inventory_name) end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = inventory%paired_object_exists(mesh_id) + constant_exists = inventory%paired_object_exists(mesh%get_id()) if (.not. constant_exists) then ! Create constant if it doesn't already exist - mass_matrix => get_mass_matrix_fv(space, mesh_id) + mass_matrix => get_mass_matrix_fv(config, mesh, space) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) @@ -871,16 +927,18 @@ contains end function get_inverse_mass_matrix_fv !> @brief Returns a pointer to the finite element curl operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @return The curl operator - function get_curl_fe(mesh_id) result(curl) + function get_curl_fe(config, mesh) result(curl) use sci_compute_curl_operator_kernel_mod, only: compute_curl_operator_kernel_type implicit none - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + type(operator_type), pointer :: curl logical(kind=l_def) :: constant_exists type(field_type), pointer :: chi(:) @@ -890,9 +948,23 @@ contains type(function_space_type), pointer :: w1_fs integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + integer(i_def) :: order_h + integer(i_def) :: order_v + real(r_def) :: scaled_radius + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + order_h = config%finite_element%element_order_h() + order_v = config%finite_element%element_order_v() + ! If running at lowest order, use finite volume - if (element_order_h == 0 .and. element_order_v == 0) then - curl => get_curl_fv(mesh_id) + if (order_h == 0 .and. order_v == 0) then + curl => get_curl_fv(config, mesh) return end if @@ -901,8 +973,7 @@ contains call curl_inventory_fe%initialise(name="curl_fe") end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = curl_inventory_fe%paired_object_exists(mesh_id) + constant_exists = curl_inventory_fe%paired_object_exists(mesh%get_id()) ! Create constant if it doesn't already exist if (.not. constant_exists) then @@ -910,16 +981,16 @@ contains panel_id => get_panel_id(mesh) qr_ptr => get_qr_fe() - w2_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W2 ) - w1_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W1 ) + w2_fs => function_space_collection%get_fs( mesh, order_h, order_v, W2 ) + w1_fs => function_space_collection%get_fs( mesh, order_h, order_v, W1 ) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) call curl_inventory_fe%add_operator(curl, w2_fs, w1_fs, mesh) - call invoke( name='calculate_curl_fe', & - compute_curl_operator_kernel_type(curl, chi, & - panel_id, qr_ptr) ) + call invoke( name='calculate_curl_fe', & + compute_curl_operator_kernel_type(curl, chi, panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, qr_ptr) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) end if @@ -928,16 +999,18 @@ contains end function get_curl_fe !> @brief Returns a pointer to the finite volume curl operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @return The curl operator - function get_curl_fv(mesh_id) result(curl) + function get_curl_fv(config, mesh) result(curl) use sci_compute_curl_operator_kernel_mod, only: compute_curl_operator_kernel_type implicit none - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + type(operator_type), pointer :: curl logical(kind=l_def) :: constant_exists type(field_type), pointer :: chi(:) @@ -947,13 +1020,22 @@ contains type(function_space_type), pointer :: w1_fs integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + ! Initialise inventory if it hasn't been done so already if (.not. curl_inventory_fv%is_initialised()) then call curl_inventory_fv%initialise(name="curl_fv") end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = curl_inventory_fv%paired_object_exists(mesh_id) + constant_exists = curl_inventory_fv%paired_object_exists(mesh%get_id()) ! Create constant if it doesn't already exist if (.not. constant_exists) then @@ -966,9 +1048,11 @@ contains if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) call curl_inventory_fv%add_operator(curl, w2_fs, w1_fs, mesh) - call invoke( name='calculate_curl_fv', & - compute_curl_operator_kernel_type(curl, chi, & - panel_id, qr_ptr) ) + call invoke( name='calculate_curl_fv', & + compute_curl_operator_kernel_type(curl, chi, panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, qr_ptr) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) end if @@ -977,14 +1061,16 @@ contains end function get_curl_fv !> @brief Returns a pointer to the div operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @return The div operator - function get_div(mesh_id) result(div) + function get_div(config, mesh) result(div) implicit none - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + type(operator_type), pointer :: div logical(kind=l_def) :: constant_exists type(field_type), pointer :: chi(:) @@ -994,13 +1080,28 @@ contains type(function_space_type), pointer :: w2_fs integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + integer(i_def) :: order_h + integer(i_def) :: order_v + real(r_def) :: scaled_radius + logical(l_def) :: rehabilitate + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + rehabilitate = config%finite_element%rehabilitate() + order_h = config%finite_element%element_order_h() + order_v = config%finite_element%element_order_v() + ! Initialise inventory if it hasn't been done so already if (.not. div_inventory%is_initialised()) then call div_inventory%initialise(name="div") end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = div_inventory%paired_object_exists(mesh_id) + constant_exists = div_inventory%paired_object_exists(mesh%get_id()) ! Create constant if it doesn't already exist if (.not. constant_exists) then @@ -1008,16 +1109,19 @@ contains panel_id => get_panel_id(mesh) qr_ptr => get_qr_fe() - w2_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W2 ) - w3_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W3 ) + w2_fs => function_space_collection%get_fs( mesh, order_h, order_v, W2 ) + w3_fs => function_space_collection%get_fs( mesh, order_h, order_v, W3 ) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) call div_inventory%add_operator(div, w3_fs, w2_fs, mesh) - call invoke( name='calculate_div', & - compute_div_operator_kernel_type(div, chi, & - panel_id, qr_ptr) ) + call invoke( name='calculate_div', & + compute_div_operator_kernel_type(div, chi, & + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + rehabilitate, & + qr_ptr) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) end if @@ -1026,14 +1130,16 @@ contains end function get_div !> @brief Returns a pointer to the div operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @return The div operator for W2H - function get_div_h(mesh_id) result(div_h) + function get_div_h(config, mesh) result(div_h) implicit none - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + type(operator_type), pointer :: div_h logical(kind=l_def) :: constant_exists type(field_type), pointer :: chi(:) @@ -1043,13 +1149,28 @@ contains type(function_space_type), pointer :: w2h_fs integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + integer(i_def) :: order_h + integer(i_def) :: order_v + real(r_def) :: scaled_radius + logical(l_def) :: rehabilitate + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + rehabilitate = config%finite_element%rehabilitate() + order_h = config%finite_element%element_order_h() + order_v = config%finite_element%element_order_v() + ! Initialise inventory if it hasn't been done so already if (.not. div_h_inventory%is_initialised()) then call div_h_inventory%initialise(name="div_h") end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = div_h_inventory%paired_object_exists(mesh_id) + constant_exists = div_h_inventory%paired_object_exists(mesh%get_id()) ! Create constant if it doesn't already exist if (.not. constant_exists) then @@ -1057,16 +1178,19 @@ contains panel_id => get_panel_id(mesh) qr_ptr => get_qr_fe() - w2h_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W2H ) - w3_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W3 ) + w2h_fs => function_space_collection%get_fs( mesh, order_h, order_v, W2H ) + w3_fs => function_space_collection%get_fs( mesh, order_h, order_v, W3 ) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) call div_h_inventory%add_operator(div_h, w3_fs, w2h_fs, mesh) - call invoke( name='calculate_div_h', & - compute_div_operator_kernel_type(div_h, chi, & - panel_id, qr_ptr) ) + call invoke( name='calculate_div_h', & + compute_div_operator_kernel_type(div_h, chi, & + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + rehabilitate, & + qr_ptr) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) end if @@ -1075,16 +1199,18 @@ contains end function get_div_h !> @brief Returns a pointer to the finite element inv_m3 * div operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @return The div operator - function get_im3_div_fe(mesh_id) result(im3_div) + function get_im3_div_fe(config, mesh) result(im3_div) use sci_operator_x_times_y_kernel_mod, only: operator_x_times_y_kernel_type implicit none - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + type(operator_type), pointer :: im3_div type(operator_type), pointer :: mm_w3_inv type(operator_type) :: div @@ -1096,9 +1222,25 @@ contains type(function_space_type), pointer :: w2_fs integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + integer(i_def) :: order_h + integer(i_def) :: order_v + real(r_def) :: scaled_radius + logical(l_def) :: rehabilitate + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + rehabilitate = config%finite_element%rehabilitate() + order_h = config%finite_element%element_order_h() + order_v = config%finite_element%element_order_v() + ! If running at lowest order, use finite volume - if (element_order_h == 0 .and. element_order_v == 0) then - im3_div => get_im3_div_fv(mesh_id) + if (order_h == 0 .and. order_v == 0) then + im3_div => get_im3_div_fv(config, mesh) return end if @@ -1107,26 +1249,28 @@ contains call im3_div_inventory_fe%initialise(name="im3_div_fe") end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = im3_div_inventory_fe%paired_object_exists(mesh_id) + constant_exists = im3_div_inventory_fe%paired_object_exists(mesh%get_id()) ! Create constant if it doesn't already exist if (.not. constant_exists) then chi => get_coordinates(mesh) panel_id => get_panel_id(mesh) qr_ptr => get_qr_fe() - mm_w3_inv => get_inverse_mass_matrix_fe(W3, mesh_id) + mm_w3_inv => get_inverse_mass_matrix_fe(config, mesh, W3) - w2_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W2 ) - w3_fs => function_space_collection%get_fs( mesh, element_order_h, & - element_order_v, W3 ) + w2_fs => function_space_collection%get_fs( mesh, order_h, order_v, W2 ) + w3_fs => function_space_collection%get_fs( mesh, order_h, order_v, W3 ) if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) call div%initialise( w3_fs, w2_fs ) call im3_div_inventory_fe%add_operator(im3_div, w3_fs, w2_fs, mesh) - call invoke( compute_div_operator_kernel_type(div, chi, & - panel_id, qr_ptr), & + call invoke( compute_div_operator_kernel_type(div, chi, & + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + rehabilitate, & + qr_ptr), & operator_x_times_y_kernel_type(im3_div, mm_w3_inv, div) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) end if @@ -1136,16 +1280,18 @@ contains end function get_im3_div_fe !> @brief Returns a pointer to the finite volume inv_m3 * div operator - !> @param[in] mesh_id The ID of the mesh to get the object for + !> @param[in] config Application configuration object + !> @param[in] mesh Mesh to get the object for !> @return The div operator - function get_im3_div_fv(mesh_id) result(im3_div) + function get_im3_div_fv(config, mesh) result(im3_div) use sci_operator_x_times_y_kernel_mod, only: operator_x_times_y_kernel_type implicit none - integer(kind=i_def), intent(in) :: mesh_id - type(mesh_type), pointer :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh + type(operator_type), pointer :: im3_div type(operator_type), pointer :: mm_w3_inv type(operator_type) :: div @@ -1157,20 +1303,31 @@ contains type(function_space_type), pointer :: w2_fs integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius + logical(l_def) :: rehabilitate + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + rehabilitate = config%finite_element%rehabilitate() + ! Initialise inventory if it hasn't been done so already if (.not. im3_div_inventory_fv%is_initialised()) then call im3_div_inventory_fv%initialise(name="im3_div_fv") end if - mesh => mesh_collection%get_mesh(mesh_id) - constant_exists = im3_div_inventory_fv%paired_object_exists(mesh_id) + constant_exists = im3_div_inventory_fv%paired_object_exists(mesh%get_id()) ! Create constant if it doesn't already exist if (.not. constant_exists) then chi => get_coordinates(mesh) panel_id => get_panel_id(mesh) qr_ptr => get_qr_fv() - mm_w3_inv => get_inverse_mass_matrix_fv(W3, mesh_id) + mm_w3_inv => get_inverse_mass_matrix_fv(config, mesh, W3) w2_fs => function_space_collection%get_fs( mesh, 0, 0, W2 ) w3_fs => function_space_collection%get_fs( mesh, 0, 0, W3 ) @@ -1178,8 +1335,13 @@ contains if ( LPROF ) call start_timing( id, 'runtime_constants.fem' ) call div%initialise( w3_fs, w2_fs ) call im3_div_inventory_fv%add_operator(im3_div, w3_fs, w2_fs, mesh) - call invoke( compute_div_operator_kernel_type(div, chi, & - panel_id, qr_ptr), & + call invoke( compute_div_operator_kernel_type(div, chi, & + panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & + rehabilitate, & + qr_ptr), & operator_x_times_y_kernel_type(im3_div, mm_w3_inv, div) ) if ( LPROF ) call stop_timing( id, 'runtime_constants.fem' ) end if diff --git a/components/science/source/algorithm/sci_galerkin_projection_alg_mod.x90 b/components/science/source/algorithm/sci_galerkin_projection_alg_mod.x90 index c726e3cf3..1ef1ea6df 100644 --- a/components/science/source/algorithm/sci_galerkin_projection_alg_mod.x90 +++ b/components/science/source/algorithm/sci_galerkin_projection_alg_mod.x90 @@ -9,6 +9,7 @@ !> module sci_galerkin_projection_alg_mod + use config_mod, only: config_type use constants_mod, only: r_def, i_def use function_space_collection_mod, only: function_space_collection use field_mod, only: field_type @@ -17,10 +18,12 @@ module sci_galerkin_projection_alg_mod use sci_gp_vector_rhs_kernel_mod, only: gp_vector_rhs_kernel_type use log_mod, only: log_event, log_scratch_space, & log_level_debug, log_level_error + use mesh_mod, only: mesh_type use operator_mod, only: operator_type use quadrature_xyoz_mod, only: quadrature_xyoz_type use sci_mass_matrix_solver_alg_mod, & only: mass_matrix_solver_alg + use sci_mesh_enums_mod, only: get_mesh_enums implicit none @@ -43,6 +46,7 @@ contains !> into orthogonal components and the galerkin projection of each !> component is computed. !> + !> @param[in] config Application configuration object !> @param[out] f_out Field to project onto. !> @param[in] f_in Field to project. !> @param[in] chi Field entity co-ordinates. @@ -52,7 +56,8 @@ contains !> @todo Ideally this function would not take chi and panel_id, these would !> be available through the input field. That is some way off though. !> - subroutine galerkin_projection_algorithm( f_out, & + subroutine galerkin_projection_algorithm( config, & + f_out, & f_in, & chi, & panel_id, & @@ -60,6 +65,8 @@ contains implicit none + type(config_type), intent(in) :: config + ! Field to output type(field_type), intent(inout) :: f_out(:) @@ -76,13 +83,27 @@ contains type(field_type) :: w2_field integer(i_def) :: out_fs integer(i_def) :: idx + + type(mesh_type), pointer :: mesh + integer(i_def) :: element_order_h integer(i_def) :: element_order_v + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius ! Create continuous fields to project data into out_fs = f_out(1)%which_function_space() + + mesh => f_in%get_mesh() element_order_h = f_in%get_element_order_h() element_order_v = f_in%get_element_order_v() + + call get_mesh_enums(mesh, geometry, topology) + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + do idx = 1, size(f_out) ! A Galerkin projection doesn't make sense between meshes. ! This check could usefully be an assertion if Fortran gains support @@ -108,7 +129,9 @@ contains if ( size(f_out) == 1 ) then write( log_scratch_space, '(A)' ) ' scalar field ... ' call log_event( log_scratch_space, log_level_debug ) - call invoke( gp_rhs_kernel_type(rhs(1), f_in, chi, panel_id, qr) ) + call invoke( gp_rhs_kernel_type(rhs(1), f_in, chi, panel_id, & + geometry, topology, coord_system, & + scaled_radius, qr) ) else write( log_scratch_space, '(A)' ) ' vector field ... ' call log_event( log_scratch_space, log_level_debug ) @@ -121,14 +144,13 @@ contains ) & ) - call invoke( & - gp_vector_rhs_kernel_type( & - rhs, f_in, chi, panel_id, w2_field, qr & - ) & - ) + call invoke( & + gp_vector_rhs_kernel_type( rhs, f_in, chi, panel_id, w2_field, & + geometry, topology, coord_system, & + scaled_radius, qr ) ) end if do idx = 1, size(f_out) - call mass_matrix_solver_alg( f_out(idx), rhs(idx) ) + call mass_matrix_solver_alg( config, f_out(idx), rhs(idx) ) end do end if diff --git a/components/science/source/algorithm/sci_geometric_constants_mod.x90 b/components/science/source/algorithm/sci_geometric_constants_mod.x90 index decdc4332..5420744d6 100644 --- a/components/science/source/algorithm/sci_geometric_constants_mod.x90 +++ b/components/science/source/algorithm/sci_geometric_constants_mod.x90 @@ -33,12 +33,9 @@ module sci_geometric_constants_mod use inventory_by_local_mesh_mod, only: inventory_by_local_mesh_type use mesh_mod, only: mesh_type use local_mesh_mod, only: local_mesh_type + use sci_mesh_enums_mod, only: get_mesh_enums, geometry_spherical ! Configuration - use base_mesh_config_mod, only: geometry_spherical, & - geometry_planar, & - topology_fully_periodic, & - topology_non_periodic use finite_element_config_mod, only: coord_system_native implicit none @@ -175,6 +172,7 @@ contains integer(kind=i_def) :: k_h, k_v integer(tik) :: id + integer(kind=i_def) :: geometry integer(kind=i_def) :: topology integer(kind=i_def) :: element_order_h integer(kind=i_def) :: element_order_v @@ -185,6 +183,8 @@ contains if ( LPROF ) call start_timing( id, 'runtime_constants.geometric' ) + call get_mesh_enums(mesh, geometry, topology) + coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() element_order_h = config%finite_element%element_order_h() @@ -214,18 +214,9 @@ contains chi => get_coordinates(mesh) panel_id => get_panel_id(mesh) - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else if (mesh%is_topology_non_periodic()) then - topology = topology_non_periodic - else - call log_event('compute_latlon: Unsupported topology selection',log_level_error) - end if - call invoke( compute_latlon_kernel_type(lat, long, chi, panel_id, & - geometry_spherical, & - topology, coord_system, & - scaled_radius) ) + geometry, topology, & + coord_system, scaled_radius) ) else call invoke( setval_c(lat, f_lat), & setval_c(long, f_lon) ) @@ -519,17 +510,27 @@ contains type(quadrature_rule_gaussian_type) :: quadrature_rule integer(tik) :: id - integer(kind=i_def) :: element_order_h, element_order_v - integer(kind=i_def) :: nqp_h_exact, nqp_v_exact - - element_order_h = config%finite_element%element_order_h() - element_order_v = config%finite_element%element_order_v() - nqp_h_exact = config%finite_element%nqp_h_exact() - nqp_v_exact = config%finite_element%nqp_v_exact() + integer(i_def) :: order_h + integer(i_def) :: order_v + integer(i_def) :: nqp_h_exact + integer(i_def) :: nqp_v_exact + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius + + call get_mesh_enums(mesh, geometry, topology) + + order_h = config%finite_element%element_order_h() + order_v = config%finite_element%element_order_v() + nqp_h_exact = config%finite_element%nqp_h_exact() + nqp_v_exact = config%finite_element%nqp_v_exact() + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() ! If running at lowest order, use finite volume - if (element_order_h == 0 .and. element_order_v == 0) then - detj_at_w3 => get_detj_at_w3_fv(mesh) + if (order_h == 0 .and. order_v == 0) then + detj_at_w3 => get_detj_at_w3_fv(config, mesh) return end if @@ -548,10 +549,7 @@ contains ! Create the object as it doesn't exist yet if ( LPROF ) call start_timing( id, 'runtime_constants.geometric' ) - w3_fs => function_space_collection%get_fs( mesh, & - element_order_h, & - element_order_v, & - W3 ) + w3_fs => function_space_collection%get_fs( mesh, order_h, order_v, W3 ) call detj_at_w3_inventory_fe%add_field(detj_at_w3, w3_fs, mesh) ! @TODO #4487: it is inefficient to calculate this via mass matrices @@ -565,6 +563,10 @@ contains chi, & panel_id, & extended_mesh, & + geometry, & + topology, & + coord_system, & + scaled_radius, & qr), & setval_c(detj_at_w3, 0.0_r_def), & mm_diagonal_kernel_type(detj_at_w3, mm_w3) ) @@ -578,9 +580,10 @@ contains end function get_detj_at_w3_fe !> @brief Returns the (finite volume) Det(J) values at W3 dof locations - !> @param[in] mesh Mesh to get the object for + !> @param[in] config Application namelist configuration object + !> @param[in] mesh Mesh to get the object for !> @return The Det(J) field - function get_detj_at_w3_fv(mesh) result(detj_at_w3) + function get_detj_at_w3_fv(config, mesh) result(detj_at_w3) ! @TODO #4487: update these imports ! use sci_calc_detj_at_w3_kernel_mod, & @@ -594,7 +597,8 @@ contains implicit none - type(mesh_type), intent(in) :: mesh + type(config_type), intent(in) :: config + type(mesh_type), intent(in) :: mesh logical(kind=l_def) :: constant_exists type(field_type), pointer :: detj_at_w3 @@ -608,6 +612,16 @@ contains type(quadrature_rule_gaussian_type) :: quadrature_rule integer(tik) :: id + integer(i_def) :: geometry + integer(i_def) :: topology + integer(i_def) :: coord_system + real(r_def) :: scaled_radius + + call get_mesh_enums(mesh, geometry, topology) + + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + ! Check inventory is initialised if (.not. detj_at_w3_inventory_fv%is_initialised()) then ! Initialise all inventories together @@ -636,6 +650,10 @@ contains chi, & panel_id, & extended_mesh, & + geometry, & + topology, & + coord_system, & + scaled_radius, & qr), & setval_c(detj_at_w3, 0.0_r_def), & mm_diagonal_kernel_type(detj_at_w3, mm_w3) ) @@ -995,11 +1013,7 @@ contains call dA_msl_proj_inventory%initialise(name="dA_msl_proj") end if - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if + call get_mesh_enums(mesh, geometry) planet_radius = config%extrusion%planet_radius() domain_height = config%extrusion%domain_height() @@ -1340,11 +1354,7 @@ contains real(kind=r_def) :: scaled_radius - if (mesh%is_geometry_spherical()) then - geometry= geometry_spherical - else - geometry= geometry_planar - end if + call get_mesh_enums(mesh, geometry) element_order_h = config%finite_element%element_order_h() element_order_v = config%finite_element%element_order_v() @@ -1470,11 +1480,7 @@ contains integer(kind=i_def) :: coord_system real(kind=r_def) :: scaled_radius - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if + call get_mesh_enums(mesh, geometry) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() diff --git a/components/science/source/algorithm/sci_mapping_constants_mod.x90 b/components/science/source/algorithm/sci_mapping_constants_mod.x90 index e295d9b75..9c1b1db4c 100644 --- a/components/science/source/algorithm/sci_mapping_constants_mod.x90 +++ b/components/science/source/algorithm/sci_mapping_constants_mod.x90 @@ -15,11 +15,11 @@ module sci_mapping_constants_mod ! Infrastructure + use config_mod, only: config_type use constants_mod, only: i_def, r_def, l_def, str_def use extrusion_mod, only: PRIME_EXTRUSION, & DOUBLE_LEVEL, & TWOD - use sci_fem_constants_mod, only: get_qr_fv use field_mod, only: field_type use function_space_mod, only: function_space_type use function_space_collection_mod, only: function_space_collection @@ -39,20 +39,10 @@ module sci_mapping_constants_mod use r_tran_field_mod, only: r_tran_field_type use timing_mod, only: start_timing, stop_timing, & tik, LPROF - - ! Object types - use config_mod, only: config_type - - ! Configuration modules - use base_mesh_config_mod, only: geometry_spherical, & - geometry_planar, & - topology_fully_periodic, & - topology_non_periodic - - - ! Other algorithms - use sci_geometric_constants_mod, only: get_coordinates, & - get_panel_id + use sci_fem_constants_mod, only: get_qr_fv + use sci_geometric_constants_mod, only: get_coordinates, & + get_panel_id + use sci_mesh_enums_mod, only: get_mesh_enums ! Kernels use sci_project_ws_to_w1_operator_kernel_mod, & @@ -230,17 +220,7 @@ contains integer(kind=i_def) :: coord_system real(kind=r_def) :: scaled_radius - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() @@ -309,17 +289,7 @@ contains integer(kind=i_def) :: coord_system real(kind=r_def) :: scaled_radius - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() @@ -696,16 +666,20 @@ contains end function get_intermesh_weights_w2 !> @brief Returns a pointer to the weights for conservative W3 mapping + !> @param[in] config Application configuration object !> @param[in] fine_mesh The fine mesh for the transform !> @param[in] coarse_mesh The coarse mesh for the transform !> @return The field containing weights for conservative W3 mapping - function get_intermesh_weights_w3_rdef(fine_mesh, coarse_mesh) result(weights_rdef) + function get_intermesh_weights_w3_rdef(config, fine_mesh, coarse_mesh) & + result(weights_rdef) use sci_fem_constants_mod, only: get_mass_matrix_diagonal_fv use sci_weights_intermesh_w3_kernel_mod, only: weights_intermesh_w3_kernel_type implicit none + type(config_type), intent(in) :: config + type(mesh_type), pointer, intent(in) :: coarse_mesh type(mesh_type), pointer, intent(in) :: fine_mesh integer(kind=i_def) :: intermesh_id @@ -741,8 +715,8 @@ contains end if ! @TODO #4487: this would be better as just the cell volumes - mm_w3_fine => get_mass_matrix_diagonal_fv(W3, fine_mesh%get_id()) - mm_w3_coarse => get_mass_matrix_diagonal_fv(W3, coarse_mesh%get_id()) + mm_w3_fine => get_mass_matrix_diagonal_fv(config, fine_mesh, W3) + mm_w3_coarse => get_mass_matrix_diagonal_fv(config, coarse_mesh, W3) if ( LPROF ) call start_timing( id, 'runtime_constants.mapping' ) @@ -768,13 +742,17 @@ contains end function get_intermesh_weights_w3_rdef !> @brief Returns a pointer to the weights for conservative W3 mapping + !> @param[in] config Application configuration object !> @param[in] fine_mesh The fine mesh for the transform !> @param[in] coarse_mesh The coarse mesh for the transform !> @return The field containing weights for conservative W3 mapping - function get_intermesh_weights_w3_rtran(fine_mesh, coarse_mesh) result(weights_rtran) + function get_intermesh_weights_w3_rtran(config, fine_mesh, coarse_mesh) & + result(weights_rtran) implicit none + type(config_type), intent(in) :: config + type(mesh_type), pointer, intent(in) :: coarse_mesh type(mesh_type), pointer, intent(in) :: fine_mesh integer(kind=i_def) :: intermesh_id @@ -800,7 +778,9 @@ contains if (.not. constant_exists) then ! Create the object as it doesn't exist yet - weights_rdef => get_intermesh_weights_w3_rdef(fine_mesh, coarse_mesh) + weights_rdef => get_intermesh_weights_w3_rdef(config, & + fine_mesh, & + coarse_mesh) if ( LPROF ) call start_timing( id, 'runtime_constants.mapping' ) @@ -1109,17 +1089,7 @@ contains integer(kind=i_def), parameter :: xdirection = 1_i_def integer(tik) :: id - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() @@ -1191,17 +1161,7 @@ contains integer(kind=i_def), parameter :: ydirection = 2_i_def integer(tik) :: id - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() @@ -1274,17 +1234,7 @@ contains integer(kind=i_def), parameter :: zdirection = 3_i_def integer(tik) :: id - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() @@ -1358,17 +1308,7 @@ contains integer(kind=i_def) :: coord_system real(kind=r_def) :: scaled_radius - if (mesh%is_geometry_spherical()) then - geometry = geometry_spherical - else - geometry = geometry_planar - end if - - if (mesh%is_topology_periodic()) then - topology = topology_fully_periodic - else - topology = topology_non_periodic - end if + call get_mesh_enums(mesh, geometry, topology) coord_system = config%finite_element%coord_system() scaled_radius = config%planet%scaled_radius() diff --git a/components/science/source/algorithm/sci_project_output_mod.F90 b/components/science/source/algorithm/sci_project_output_mod.F90 index 8074c742d..6783b6409 100644 --- a/components/science/source/algorithm/sci_project_output_mod.F90 +++ b/components/science/source/algorithm/sci_project_output_mod.F90 @@ -20,16 +20,18 @@ module sci_project_output_mod !> @details This procedure uses the galerkin projection and a precomputed !> mass matrix to project a field !> + !> @param[in] config Application configuration object !> @param[in] field To be projected. !> @param[inout] projected_field Receives projection. !> @param[in] chi Field entity co-ordinates. !> @param[in] panel_id Cell orientation map. !> @param[in] output_fs Desired output function space. !> - subroutine project_output( field, projected_field, & + subroutine project_output( config, field, projected_field, & chi, panel_id, & output_fs ) + use config_mod, only: config_type use constants_mod, only: r_def, str_max_filename, i_def use field_mod, only: field_type use field_parent_mod, only: write_interface @@ -42,13 +44,18 @@ subroutine project_output( field, projected_field, & implicit none + type(config_type), intent(in) :: config + ! Input field to project from type(field_type), intent(in) :: field + ! Output field to project to type(field_type), intent(inout) :: projected_field(:) + ! Co-ordinate system type(field_type), intent(in) :: chi(:) type(field_type), intent(in) :: panel_id + ! Output function space integer(i_def), intent(in) :: output_fs @@ -84,9 +91,8 @@ subroutine project_output( field, projected_field, & end do ! do the projection - call galerkin_projection_algorithm( & - projected_field, field, chi, panel_id, qr & - ) + call galerkin_projection_algorithm(config, projected_field, field, & + chi, panel_id, qr) end subroutine project_output diff --git a/components/science/source/algorithm/solver/sci_dense_operator_alg_mod.x90 b/components/science/source/algorithm/solver/sci_dense_operator_alg_mod.x90 index de9e14b18..5b34d14ab 100644 --- a/components/science/source/algorithm/solver/sci_dense_operator_alg_mod.x90 +++ b/components/science/source/algorithm/solver/sci_dense_operator_alg_mod.x90 @@ -21,7 +21,8 @@ module sci_dense_operator_alg_mod only : compute_mass_matrix_kernel_w_scalar_type use sci_compute_mass_matrix_kernel_w3_mod, & only : compute_mass_matrix_kernel_w3_type - use constants_mod, only : r_def, i_def, l_def + use config_mod, only : config_type + use constants_mod, only : r_def, i_def, l_def, imdi use field_mod, only : field_type use function_space_mod, only : function_space_type use function_space_collection_mod, only : function_space_collection @@ -36,6 +37,7 @@ module sci_dense_operator_alg_mod use quadrature_xyoz_mod, only : quadrature_xyoz_type use sci_field_vector_mod, only : field_vector_type use sci_linear_operator_mod, only : abstract_linear_operator_type + use sci_mesh_enums_mod, only : get_mesh_enums use vector_mod, only : abstract_vector_type implicit none @@ -71,15 +73,18 @@ module sci_dense_operator_alg_mod contains !> @brief Construct a dense_operator object. !> + !> @param[in] config Application configuration object !> @param[in] fv A field_vector of size 2 !> @param[in] chi The coordinate field !> @param[in] panel_id The field with mesh panel IDs !> @param[in] qr The quadrature object !> @return self The constructed dense operator - function dense_operator_constructor( fv, qr, chi, panel_id ) result(self) + function dense_operator_constructor( config, fv, qr, chi, panel_id ) result(self) implicit none + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(in) :: fv type( field_type ), dimension(3), intent(in) :: chi type( field_type ), intent(in) :: panel_id @@ -87,13 +92,23 @@ contains type(dense_operator_type) :: self + type(mesh_type), pointer :: mesh + type(function_space_type), pointer :: fs_wa => null() type(function_space_type), pointer :: fs_wb => null() - type(mesh_type), pointer :: mesh => null() + integer(kind=i_def) :: element_order_h, element_order_v, fs_label logical(kind=l_def) :: extended_mesh + logical(l_def) :: rehabilitate + integer(i_def) :: geometry, topology, coord_system + real(r_def) :: scaled_radius + + rehabilitate = config%finite_element%rehabilitate() + coord_system = config%finite_element%coord_system() + scaled_radius = config%planet%scaled_radius() + select type (fv) type is (field_vector_type) @@ -105,6 +120,7 @@ contains end if mesh => fv%vector(1)%get_mesh() + call get_mesh_enums(mesh, geometry, topology) ! Make the mass_w0 operator element_order_h = fv%vector(1)%get_element_order_h() @@ -125,8 +141,12 @@ contains extended_mesh = .false. call invoke( & compute_mass_matrix_kernel_w_scalar_type(self%mass_w0, chi, panel_id, & - extended_mesh, qr), & - compute_mass_matrix_kernel_w3_type(self%mass_w3, chi, panel_id, qr) ) + extended_mesh, geometry, & + topology, coord_system, & + scaled_radius, qr), & + compute_mass_matrix_kernel_w3_type(self%mass_w3, chi, panel_id, & + geometry, topology, coord_system, & + scaled_radius, rehabilitate, qr) ) class default @@ -141,15 +161,19 @@ contains !> @brief Applies the LMA mass matrices to the vector, y = M x. !> - !> @param[in,out] self The dense operator - !> @param[in] x Field vector to be read, size 2 - !> @param[in,out] y Field vector to be written, size 2 - subroutine apply_dense_op(self, x, y) + !> @param[in] config Application configuration object + !> @param[in] x Field vector to be read, size 2 + !> @param[in,out] y Field vector to be written, size 2 + subroutine apply_dense_op(self, config, x, y) implicit none - class(dense_operator_type), intent(inout) :: self - class(abstract_vector_type), intent(in) :: x - class(abstract_vector_type), intent(inout) :: y + + class(dense_operator_type), intent(inout) :: self + + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(in) :: x + class(abstract_vector_type), intent(inout) :: y + type(field_type), pointer :: x_field1 => null(), & y_field1 => null(), & x_field2 => null(), & diff --git a/components/science/source/algorithm/solver/sci_diagonal_preconditioner_alg_mod.x90 b/components/science/source/algorithm/solver/sci_diagonal_preconditioner_alg_mod.x90 index a04d49a95..fcaad004e 100644 --- a/components/science/source/algorithm/solver/sci_diagonal_preconditioner_alg_mod.x90 +++ b/components/science/source/algorithm/solver/sci_diagonal_preconditioner_alg_mod.x90 @@ -7,6 +7,7 @@ !> @brief Application of the diagonal of an operator. module sci_diagonal_preconditioner_alg_mod + use config_mod, only: config_type use constants_mod, only: i_def, r_def use field_mod, only: field_type use log_mod, only: log_event, & @@ -65,16 +66,18 @@ contains !> @brief Apply diagonal preconditioner to a field to obtain \f$y=D^-1x\f$. !> - !> @param[in] self Instance of diagonal_preconditioner_type - !> @param[in] x Field \f$x\f$ to apply preconditioner to - !> @param[in,out] y Resulting field \f$y=D^-1x\f$ - subroutine apply_diagonal_preconditioner(self, x, y) + !> @param[in] config Application configuration object + !> @param[in] x Field \f$x\f$ to apply preconditioner to + !> @param[in,out] y Resulting field \f$y=D^-1x\f$ + subroutine apply_diagonal_preconditioner(self, config, x, y) implicit none class(diagonal_preconditioner_type), intent(inout) :: self - class(abstract_vector_type), intent(in) :: x - class(abstract_vector_type), intent(inout) :: y + + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(in) :: x + class(abstract_vector_type), intent(inout) :: y integer(kind=i_def) :: i, nfields ! Workaround for PSyclone to get pointers of the correct type for x and y diff --git a/components/science/source/algorithm/solver/sci_hori_mass_matrix_solver_alg_mod.x90 b/components/science/source/algorithm/solver/sci_hori_mass_matrix_solver_alg_mod.x90 index 366d22313..a5a4f31b8 100644 --- a/components/science/source/algorithm/solver/sci_hori_mass_matrix_solver_alg_mod.x90 +++ b/components/science/source/algorithm/solver/sci_hori_mass_matrix_solver_alg_mod.x90 @@ -12,6 +12,7 @@ module sci_hori_mass_matrix_solver_alg_mod ! Constants and types + use config_mod, only: config_type use constants_mod, only: i_def use integer_field_mod, only: integer_field_type use field_mod, only: field_type @@ -35,6 +36,7 @@ contains ! ======================================================= !> @brief Mass matrix solver for only the horizontal component of a W2 field + !> @param[in] config Application configuration object !> @param[in,out] field_solve The W2 solved field !> @param[in] field The W2 field to apply the horizontal mass !! matrix solve @@ -42,10 +44,12 @@ contains !! not inputted then the vertical solved !! component is set to zero - subroutine hori_mass_matrix_solver_alg( field_solve, field, vert_field_solve ) + subroutine hori_mass_matrix_solver_alg( config, field_solve, field, vert_field_solve ) implicit none + type(config_type), intent(in) :: config + ! Arguments type(field_type), intent(inout) :: field_solve type(field_type), intent(in) :: field @@ -76,7 +80,7 @@ contains ! Take horizontal part of field call split_w2_field_alg( hori_field, vert_field, field ) ! Mass matrix solve of horizontal part - call mass_matrix_solver_alg(hori_field_solve, hori_field) + call mass_matrix_solver_alg(config, hori_field_solve, hori_field) ! Convert to W2 field if (present(vert_field_solve)) then @@ -100,4 +104,4 @@ contains end subroutine hori_mass_matrix_solver_alg -end module sci_hori_mass_matrix_solver_alg_mod \ No newline at end of file +end module sci_hori_mass_matrix_solver_alg_mod diff --git a/components/science/source/algorithm/solver/sci_mass_matrix_operator_alg_mod.x90 b/components/science/source/algorithm/solver/sci_mass_matrix_operator_alg_mod.x90 index 8a56ff1e0..928ea0854 100644 --- a/components/science/source/algorithm/solver/sci_mass_matrix_operator_alg_mod.x90 +++ b/components/science/source/algorithm/solver/sci_mass_matrix_operator_alg_mod.x90 @@ -9,6 +9,7 @@ !> @details Applies the appropriate mass matrix operator to a field. module sci_mass_matrix_operator_alg_mod + use config_mod, only: config_type use constants_mod, only: i_def, r_def, l_def use field_mod, only: field_type use log_mod, only: log_event, & @@ -66,10 +67,10 @@ contains !> @brief Apply mass matrix operator to a pressure field to obtain \f$y=Mx\f$. !> - !> @param[in,out] self Instance of mass_matrix_operator_type - !> @param[in] x Field \f$x\f$ to apply operator to - !> @param[in,out] y Resulting field \f$y=Mx\f$ - subroutine apply_mass_matrix_operator(self, x, y) + !> @param[in] config Application configuration object + !> @param[in] x Field \f$x\f$ to apply operator to + !> @param[in,out] y Resulting field \f$y=Mx\f$ + subroutine apply_mass_matrix_operator(self, config, x, y) use sci_enforce_bc_kernel_mod, only: enforce_bc_kernel_type use matrix_vector_kernel_mod, only: matrix_vector_kernel_type @@ -77,8 +78,10 @@ contains implicit none class(mass_matrix_operator_type), intent(inout) :: self - class(abstract_vector_type), intent(in) :: x - class(abstract_vector_type), intent(inout) :: y + + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(in) :: x + class(abstract_vector_type), intent(inout) :: y ! Workaround for PSyclone to get pointers of the correct type for x and y type(field_type), pointer :: x_vec => null(), & diff --git a/components/science/source/algorithm/solver/sci_mass_matrix_solver_alg_mod.x90 b/components/science/source/algorithm/solver/sci_mass_matrix_solver_alg_mod.x90 index 20c1246b8..1c1ce9a4d 100644 --- a/components/science/source/algorithm/solver/sci_mass_matrix_solver_alg_mod.x90 +++ b/components/science/source/algorithm/solver/sci_mass_matrix_solver_alg_mod.x90 @@ -14,12 +14,14 @@ !> created when needed module sci_mass_matrix_solver_alg_mod + use config_mod, only: config_type use constants_mod, only: i_def, r_def, l_def use log_mod, only: log_event, & log_scratch_space, & LOG_LEVEL_DEBUG, & LOG_LEVEL_ERROR, & LOG_LEVEL_TRACE + use mesh_mod, only: mesh_type ! Derived Types use field_mod, only: field_type @@ -54,11 +56,13 @@ module sci_mass_matrix_solver_alg_mod contains !=============================================================================! !> @details Combined initialisation, run and finalise procedure for the mass matrix solver to solve My = x + !> @param[in] config Application configuration object !> @param[inout] y result field !> @param[in] x input field !> @param[in] bc_flag optional flag to overide the default application of !> boundary conditions - subroutine mass_matrix_solver_alg(y, x, bc_flag) + subroutine mass_matrix_solver_alg(config, y, x, bc_flag) + use solver_config_mod, only: maximum_iterations, & tolerance, & method, & @@ -91,9 +95,12 @@ contains only: operator_tri_solve_kernel_type implicit none + type(config_type), intent(in) :: config + ! Prognostic fields - type( field_type), intent(inout) :: y - type( field_type), intent(in) :: x + type( field_type), intent(inout) :: y + type( field_type), intent(in) :: x + ! Optional flag to overide default boundary condition application logical(kind=l_def), optional, intent(in) :: bc_flag @@ -111,7 +118,9 @@ contains integer(kind=i_def) :: fs logical(kind=l_def) :: apply_bc logical(kind=l_def) :: lowest_order - integer(kind=i_def) :: mesh_id + + type(mesh_type), pointer :: mesh + type(operator_type), pointer :: mass_matrix => null() type(field_type), pointer :: mass_matrix_diagonal => null() type(field_vector_type) :: vec_mm_diagonal @@ -121,7 +130,7 @@ contains if ( LPROF ) call start_timing( id, 'mass_matrix_solver_alg' ) - mesh_id = y%get_mesh_id() + mesh => y%get_mesh() lowest_order = x%get_element_order_h() == 0 .and. & x%get_element_order_v() == 0 @@ -131,14 +140,14 @@ contains ! Direct solve for fully discontinuous spaces !> todo this needs to be extended to all DG spaces if ( lowest_order ) then - mass_matrix => get_inverse_mass_matrix_fv(fs, mesh_id) + mass_matrix => get_inverse_mass_matrix_fv(config, mesh, fs) else - mass_matrix => get_inverse_mass_matrix_fe(fs, mesh_id) + mass_matrix => get_inverse_mass_matrix_fe(config, mesh, fs) end if call invoke( dg_matrix_vector_kernel_type(y, x, mass_matrix) ) else if ( fs == Wtheta .and. lowest_order ) then - mass_matrix => get_mass_matrix_fv(fs, mesh_id) + mass_matrix => get_mass_matrix_fv(config, mesh, fs) call invoke( operator_tri_solve_kernel_type(y, x, mass_matrix) ) else ! Iterative solve for (semi-)continuous spaces @@ -192,12 +201,12 @@ contains if ( lowest_order ) then ! Use mass matrix operator for FV spaces - mass_matrix => get_mass_matrix_fv(fs, mesh_id) - mass_matrix_diagonal => get_mass_matrix_diagonal_fv(fs, mesh_id) + mass_matrix => get_mass_matrix_fv(config, mesh, fs) + mass_matrix_diagonal => get_mass_matrix_diagonal_fv(config, mesh, fs) else ! Use mass matrix operator for higher order spaces - mass_matrix => get_mass_matrix_fe(fs, mesh_id) - mass_matrix_diagonal => get_mass_matrix_diagonal_fe(fs, mesh_id) + mass_matrix => get_mass_matrix_fe(config, mesh, fs) + mass_matrix_diagonal => get_mass_matrix_diagonal_fe(config, mesh, fs) end if select case ( fs ) @@ -307,7 +316,7 @@ contains vec_y = field_vector_type(1) call vec_x%import_field(x, 1) call vec_y%import_field(y, 1) - call mass_matrix_solver%apply(vec_y, vec_x) + call mass_matrix_solver%apply(config, vec_y, vec_x) call vec_y%export_field(y, 1) end if ! direct or iterative solve diff --git a/components/science/source/algorithm/solver/sci_null_preconditioner_alg_mod.x90 b/components/science/source/algorithm/solver/sci_null_preconditioner_alg_mod.x90 index c229fff55..e2e9b6321 100644 --- a/components/science/source/algorithm/solver/sci_null_preconditioner_alg_mod.x90 +++ b/components/science/source/algorithm/solver/sci_null_preconditioner_alg_mod.x90 @@ -11,6 +11,7 @@ module sci_null_preconditioner_alg_mod + use config_mod, only : config_type use constants_mod, only : i_def, r_def use log_mod, only : log_event, & LOG_LEVEL_INFO, & @@ -67,15 +68,18 @@ contains !> !>@details Apply the trivial (null) preconditioner by setting \f$y=x\f$ !> - !>@param[inout] self instance of type null_preconditioner_type - !>@param[inout] x field-vector containing the right hand side of the pressure - !>@param[inout] y field-vector containing the solution - subroutine apply_null_preconditioner(self, x, y) + !>@param[in] config Application configuration object + !>@param[inout] x field-vector containing the right hand side of the pressure + !>@param[inout] y field-vector containing the solution + subroutine apply_null_preconditioner(self, config, x, y) implicit none + class(null_preconditioner_type), intent(inout) :: self - class(abstract_vector_type), intent(in) :: x - class(abstract_vector_type), intent(inout) :: y + + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(in) :: x + class(abstract_vector_type), intent(inout) :: y select type(x) type is(field_vector_type) diff --git a/components/science/source/kernel/fem/sci_compute_broken_div_operator_kernel_mod.F90 b/components/science/source/kernel/fem/sci_compute_broken_div_operator_kernel_mod.F90 index 3d00ff143..f09402ef7 100644 --- a/components/science/source/kernel/fem/sci_compute_broken_div_operator_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_broken_div_operator_kernel_mod.F90 @@ -7,21 +7,19 @@ module sci_compute_broken_div_operator_kernel_mod use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_1, & + GH_REAL, GH_INTEGER, & + GH_LOGICAL, ANY_SPACE_1, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ - use constants_mod, only: r_def, i_def + use constants_mod, only: r_def, i_def, l_def use sci_coordinate_jacobian_mod, only: coordinate_jacobian use fs_continuity_mod, only: W2broken, W3 use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system, rehabilitate - use planet_config_mod, only: scaled_radius - implicit none private @@ -32,10 +30,15 @@ module sci_compute_broken_div_operator_kernel_mod type, public, extends(kernel_type) :: compute_broken_div_operator_kernel_type private - type(arg_type) :: meta_args(3) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2broken), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(8) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2broken), & ! broken_div + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ), & ! scaled_radius + arg_type(GH_SCALAR, GH_LOGICAL, GH_READ) & ! rehabilitate /) type(func_type) :: meta_funcs(3) = (/ & func_type(W3, GH_BASIS), & @@ -64,6 +67,11 @@ module sci_compute_broken_div_operator_kernel_mod !! @param[in] chi2 2nd coordinate field in Wchi !! @param[in] chi3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels. + !! @param[in] geometry Mesh geometry enumeration + !! @param[in] topology Mesh topology enumeration + !! @param[in] coord_system Finite-element coordinate system enumeration + !! @param[in] scaled_radius Scaled planet radius + !! @param[in] rehabilitate Apply rehabilitation !! @param[in] ndf_w3 Number of degrees of freedom per cell for W3 space. !! @param[in] basis_w3 Scalar basis functions evaluated at quadrature points !! for W3 space. @@ -89,6 +97,9 @@ module sci_compute_broken_div_operator_kernel_mod subroutine compute_broken_div_operator_code(cell, nlayers, ncell_3d, & broken_div, & chi1, chi2, chi3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & + rehabilitate, & ndf_w3, basis_w3, & ndf_w2b, diff_basis_w2b, & ndf_chi, undf_chi, map_chi, & @@ -122,6 +133,12 @@ subroutine compute_broken_div_operator_code(cell, nlayers, ncell_3d, & real(kind=r_def), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_def), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + logical(kind=l_def), intent(in) :: rehabilitate + ! Internal variables integer(kind=i_def) :: df, df2, df3, k, ik, ipanel integer(kind=i_def) :: qp1, qp2 diff --git a/components/science/source/kernel/fem/sci_compute_curl_operator_kernel_mod.F90 b/components/science/source/kernel/fem/sci_compute_curl_operator_kernel_mod.F90 index ddcb67d05..b504249fd 100644 --- a/components/science/source/kernel/fem/sci_compute_curl_operator_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_curl_operator_kernel_mod.F90 @@ -6,21 +6,20 @@ module sci_compute_curl_operator_kernel_mod use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_1, & + GH_REAL, GH_INTEGER, & + ANY_SPACE_1, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ + use constants_mod, only: r_def, i_def use sci_coordinate_jacobian_mod, only: coordinate_jacobian use fs_continuity_mod, only: W1, W2 use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none private @@ -31,10 +30,14 @@ module sci_compute_curl_operator_kernel_mod type, public, extends(kernel_type) :: compute_curl_operator_kernel_type private - type(arg_type) :: meta_args(3) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W1), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(7) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W1), & ! curl + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(3) = (/ & func_type(W2, GH_BASIS), & @@ -63,6 +66,10 @@ module sci_compute_curl_operator_kernel_mod !! @param[in] chi2 2nd component of coordinate field !! @param[in] chi3 3rd component of coordinate field !! @param[in] panel_id A field giving the ID for mesh panels. +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf_w2 Number of degrees of freedom per cell for W2. !! @param[in] basis_w2 W2 vector basis functions evaluated at quadrature points. !! @param[in] ndf_w1 Number of degrees of freedom per cell for W1. @@ -83,6 +90,8 @@ module sci_compute_curl_operator_kernel_mod subroutine compute_curl_operator_code(cell, nlayers, ncell_3d, & curl, & chi1, chi2, chi3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w2, basis_w2, & ndf_w1, diff_basis_w1, & ndf_chi, undf_chi, map_chi, & @@ -115,6 +124,11 @@ subroutine compute_curl_operator_code(cell, nlayers, ncell_3d, & real(kind=r_def), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_def), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df1, df2, k, ik, ipanel integer(kind=i_def) :: qp1, qp2 diff --git a/components/science/source/kernel/fem/sci_compute_derham_matrices_kernel_mod.F90 b/components/science/source/kernel/fem/sci_compute_derham_matrices_kernel_mod.F90 index dcbf2607b..45ada61c4 100644 --- a/components/science/source/kernel/fem/sci_compute_derham_matrices_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_derham_matrices_kernel_mod.F90 @@ -16,9 +16,11 @@ module sci_compute_derham_matrices_kernel_mod use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & - GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_9, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & + GH_READ, GH_WRITE, & + GH_REAL, GH_INTEGER, & + ANY_SPACE_9, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ @@ -28,10 +30,6 @@ module sci_compute_derham_matrices_kernel_mod use fs_continuity_mod, only: W0, W1, W2, W2broken, W3, Wtheta use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none private @@ -42,19 +40,23 @@ module sci_compute_derham_matrices_kernel_mod type, public, extends(kernel_type) :: compute_derham_matrices_kernel_type private - type(arg_type) :: meta_args(12) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W0, W0), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W1), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W2), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2broken, W2broken), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W3), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, Wtheta, Wtheta), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W0), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W1), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2broken), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(16) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W0, W0), & ! mm0 + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W1), & ! mm1 + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W2), & ! mm2 + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2broken, W2broken), & ! mm2b + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W3), & ! mm3 + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, Wtheta, Wtheta), & ! mmt + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W0), & ! grad + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W1), & ! curl + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2), & ! div + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2broken), & ! broken_div + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(7) = (/ & func_type(W0, GH_BASIS, GH_DIFF_BASIS), & @@ -105,6 +107,11 @@ module sci_compute_derham_matrices_kernel_mod !! @param[in] chi2 Physical coordinates in the 2nd dir. !! @param[in] chi3 Physical coordinates in the 3rd dir. !! @param[in] panel_id Field giving the ID for mesh panels. +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius +!! @param[in] rehabilitate Apply rehabilitation !! @param[in] ndf_w0 Number of degrees of freedom per cell for W0 space. !! @param[in] basis_w0 Basis functions evaluated at quadrature points for W0 space. !! @param[in] diff_basis_w0 Differential of basis functions evaluated at quadrature points for W0 space. @@ -148,6 +155,8 @@ subroutine compute_derham_matrices_code(cell, nlayers, & ncell_3d7, broken_div, & chi1, chi2, chi3, & panel_id, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w0, basis_w0, diff_basis_w0, & ndf_w1, basis_w1, diff_basis_w1, & ndf_w2, basis_w2, diff_basis_w2, & @@ -205,6 +214,11 @@ subroutine compute_derham_matrices_code(cell, nlayers, & real(kind=r_def), intent(in) :: wqp_h(nqp_h) real(kind=r_def), intent(in) :: wqp_v(nqp_v) + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df2, k, ik integer(kind=i_def) :: qp1, qp2 diff --git a/components/science/source/kernel/fem/sci_compute_div_operator_kernel_mod.F90 b/components/science/source/kernel/fem/sci_compute_div_operator_kernel_mod.F90 index 05b84b344..93eec65a7 100644 --- a/components/science/source/kernel/fem/sci_compute_div_operator_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_div_operator_kernel_mod.F90 @@ -6,21 +6,19 @@ module sci_compute_div_operator_kernel_mod use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_1, & + GH_REAL, GH_INTEGER, & + GH_LOGICAL, ANY_SPACE_1, & GH_BASIS, GH_DIFF_BASIS, & ANY_DISCONTINUOUS_SPACE_3, & GH_QUADRATURE_XYoZ, CELL_COLUMN - use constants_mod, only: r_def, i_def + use constants_mod, only: r_def, i_def, l_def use sci_coordinate_jacobian_mod, only: coordinate_jacobian use fs_continuity_mod, only: W2, W3 use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system, rehabilitate - use planet_config_mod, only: scaled_radius - implicit none private @@ -31,10 +29,16 @@ module sci_compute_div_operator_kernel_mod type, public, extends(kernel_type) :: compute_div_operator_kernel_type private - type(arg_type) :: meta_args(3) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(8) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2), & ! div + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ), & ! scaled_radius + arg_type(GH_SCALAR, GH_LOGICAL, GH_READ) & ! rehabilitate + /) type(func_type) :: meta_funcs(3) = (/ & func_type(W3, GH_BASIS), & @@ -63,6 +67,11 @@ module sci_compute_div_operator_kernel_mod !! @param[in] chi_2 2nd coordinate field in Wchi !! @param[in] chi_3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels. +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius +!! @param[in] rehabilitate Apply rehabilitation !! @param[in] ndf_w3 Number of degrees of freedom per cell. !! @param[in] basis_w3 Scalar basis functions !! evaluated at quadrature points. @@ -87,6 +96,9 @@ subroutine compute_div_operator_code(cell, nlayers, ncell_3d, & div, & chi1, chi2, chi3, & panel_id, & + geometry, topology, & + coord_system, scaled_radius, & + rehabilitate, & ndf_w3, & basis_w3, & ndf_w2, diff_basis_w2, & @@ -124,6 +136,12 @@ subroutine compute_div_operator_code(cell, nlayers, ncell_3d, & real(kind=r_def), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_def), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + logical(kind=l_def), intent(in) :: rehabilitate + ! Internal variables integer(kind=i_def) :: df, df2, df3, k, ik integer(kind=i_def) :: qp1, qp2, ipanel diff --git a/components/science/source/kernel/fem/sci_compute_grad_operator_kernel_mod.F90 b/components/science/source/kernel/fem/sci_compute_grad_operator_kernel_mod.F90 index 9c06d929f..2006ebc4a 100644 --- a/components/science/source/kernel/fem/sci_compute_grad_operator_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_grad_operator_kernel_mod.F90 @@ -6,9 +6,11 @@ module sci_compute_grad_operator_kernel_mod use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_1, & + GH_REAL, GH_INTEGER, & + ANY_SPACE_1, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ @@ -19,10 +21,6 @@ module sci_compute_grad_operator_kernel_mod use fs_continuity_mod, only: W0, W1 use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none private @@ -33,10 +31,14 @@ module sci_compute_grad_operator_kernel_mod type, public, extends(kernel_type) :: compute_grad_operator_kernel_type private - type(arg_type) :: meta_args(3) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W0), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(7) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W0), & ! grad + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_1), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(3) = (/ & func_type(W1, GH_BASIS), & @@ -65,6 +67,10 @@ module sci_compute_grad_operator_kernel_mod !! @param[in] chi2 2nd coordinate field in Wchi !! @param[in] chi3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf_w1 Number of degrees of freedom per cell !! @param[in] basis_w1 Vector basis functions !! evaluated at quadrature points @@ -85,14 +91,16 @@ module sci_compute_grad_operator_kernel_mod !! @param[in] nqp_v Number of vertical quadrature points !! @param[in] wqp_h Horizontal quadrature weights !! @param[in] wqp_v Vertical quadrature weights -subroutine compute_grad_operator_code(cell, nlayers, ncell_3d, & - grad, & - chi1, chi2, chi3, panel_id, & - ndf_w1, basis_w1, & - ndf_w0, diff_basis_w0, & - ndf_chi, undf_chi, map_chi, & - basis_chi, diff_basis_chi, & - ndf_pid, undf_pid, map_pid, & +subroutine compute_grad_operator_code(cell, nlayers, ncell_3d, & + grad, & + chi1, chi2, chi3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & + ndf_w1, basis_w1, & + ndf_w0, diff_basis_w0, & + ndf_chi, undf_chi, map_chi, & + basis_chi, diff_basis_chi, & + ndf_pid, undf_pid, map_pid, & nqp_h, nqp_v, wqp_h, wqp_v ) implicit none @@ -120,6 +128,11 @@ subroutine compute_grad_operator_code(cell, nlayers, ncell_3d, & real(kind=r_def), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_def), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df0, df1, k, ik integer(kind=i_def) :: qp1, qp2, ipanel diff --git a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w1_mod.F90 b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w1_mod.F90 index b0beabbbb..045c228fd 100644 --- a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w1_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w1_mod.F90 @@ -11,22 +11,21 @@ module sci_compute_mass_matrix_kernel_w1_mod use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_9, & + GH_REAL, GH_INTEGER, & + ANY_SPACE_9, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ + use sci_coordinate_jacobian_mod, only: coordinate_jacobian, & coordinate_jacobian_inverse use constants_mod, only: r_def, i_def use fs_continuity_mod, only: W1 use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none private @@ -37,10 +36,14 @@ module sci_compute_mass_matrix_kernel_w1_mod type, public, extends(kernel_type) :: compute_mass_matrix_kernel_w1_type private - type(arg_type) :: meta_args(3) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W1), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(7) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W1, W1), & ! mm + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(2) = (/ & func_type(ANY_SPACE_9, GH_BASIS, GH_DIFF_BASIS), & @@ -68,6 +71,10 @@ module sci_compute_mass_matrix_kernel_w1_mod !! @param[in] chi2 2nd coordinate field in Wchi !! @param[in] chi3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf_w1 Number of degrees of freedom per cell !! @param[in] basis_w1 Vector basis functions evaluated at quadrature points !! @param[in] ndf_chi Number of degrees of freedom per cell for chi field @@ -88,6 +95,8 @@ module sci_compute_mass_matrix_kernel_w1_mod subroutine compute_mass_matrix_w1_code(cell, nlayers, ncell_3d, & mm, & chi1, chi2, chi3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w1, basis_w1, & ndf_chi, undf_chi, map_chi, & basis_chi, diff_basis_chi, & @@ -116,6 +125,11 @@ subroutine compute_mass_matrix_w1_code(cell, nlayers, ncell_3d, & real(kind=r_def), intent(in) :: wqp_h(nqp_h) real(kind=r_def), intent(in) :: wqp_v(nqp_v) + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df2, k, ik, ipanel integer(kind=i_def) :: qp1, qp2 diff --git a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w2_mod.F90 b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w2_mod.F90 index f1ebbcaae..a97369886 100644 --- a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w2_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w2_mod.F90 @@ -10,22 +10,19 @@ !> module sci_compute_mass_matrix_kernel_w2_mod - use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & - GH_READ, GH_WRITE, & - GH_REAL, ANY_W2, & - ANY_DISCONTINUOUS_SPACE_3, & - ANY_SPACE_9, & - GH_BASIS, GH_DIFF_BASIS, & + use argument_mod, only: arg_type, func_type, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & + GH_READ, GH_WRITE, & + GH_REAL, GH_INTEGER, ANY_W2, & + ANY_DISCONTINUOUS_SPACE_3, & + ANY_SPACE_9, & + GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ use constants_mod, only: i_def, r_def use sci_coordinate_jacobian_mod, only: coordinate_jacobian use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none private @@ -36,10 +33,14 @@ module sci_compute_mass_matrix_kernel_w2_mod type, public, extends(kernel_type) :: compute_mass_matrix_kernel_w2_type private - type(arg_type) :: meta_args(3) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, ANY_W2, ANY_W2), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(7) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, ANY_W2, ANY_W2), & ! mm + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(2) = (/ & func_type(ANY_W2, GH_BASIS), & @@ -68,6 +69,10 @@ module sci_compute_mass_matrix_kernel_w2_mod !! @param[in] chi2 2nd coordinate field in Wchi !! @param[in] chi3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf_w2 Degrees of freedom per cell !! @param[in] basis_w2 Vector basis functions evaluated at quadrature points !! @param[in] ndf_chi Degrees of freedom per cell for chi field @@ -88,6 +93,8 @@ subroutine compute_mass_matrix_w2_code(cell, nlayers, ncell_3d, & mm, & chi1, chi2, chi3, & panel_id, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w2, basis_w2, & ndf_chi, undf_chi, map_chi, & basis_chi, & @@ -120,6 +127,11 @@ subroutine compute_mass_matrix_w2_code(cell, nlayers, ncell_3d, & real(kind=r_def), intent(in) :: wqp_h(nqp_h) real(kind=r_def), intent(in) :: wqp_v(nqp_v) + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df2, k, ik, ipanel integer(kind=i_def) :: qp1, qp2 diff --git a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w3_mod.F90 b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w3_mod.F90 index ac030d6be..27590f8bc 100644 --- a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w3_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w3_mod.F90 @@ -9,21 +9,19 @@ !> module sci_compute_mass_matrix_kernel_w3_mod - use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & - GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_9, & - ANY_DISCONTINUOUS_SPACE_3, & - GH_BASIS, GH_DIFF_BASIS, & - CELL_COLUMN, GH_QUADRATURE_XYoZ + use argument_mod, only: arg_type, func_type, & + GH_OPERATOR, GH_FIELD, GH_SCALAR, & + GH_READ, GH_WRITE, & + GH_REAL, GH_INTEGER, GH_LOGICAL, & + ANY_SPACE_9, & + ANY_DISCONTINUOUS_SPACE_3, & + GH_BASIS, GH_DIFF_BASIS, & + CELL_COLUMN, GH_QUADRATURE_XYoZ use sci_coordinate_jacobian_mod, only: coordinate_jacobian - use constants_mod, only: i_def, r_single, r_double - use fs_continuity_mod, only: W3 - use kernel_mod, only: kernel_type - - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system, rehabilitate - use planet_config_mod, only: scaled_radius + use constants_mod, only: i_def, r_def, l_def, & + r_single, r_double + use fs_continuity_mod, only: W3 + use kernel_mod, only: kernel_type implicit none @@ -34,10 +32,15 @@ module sci_compute_mass_matrix_kernel_w3_mod !--------------------------------------------------------------------------- type, public, extends(kernel_type) :: compute_mass_matrix_kernel_w3_type private - type(arg_type) :: meta_args(3) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W3), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(8) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W3), & ! mm + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ), & ! scaled_radius + arg_type(GH_SCALAR, GH_LOGICAL, GH_READ) & ! rehabilitate /) type(func_type) :: meta_funcs(2) = (/ & func_type(W3, GH_BASIS), & @@ -70,6 +73,11 @@ module sci_compute_mass_matrix_kernel_w3_mod !! @param[in] chi2 2nd coordinate field in Wchi !! @param[in] chi3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels + !! @param[in] geometry Mesh geometry enumeration + !! @param[in] topology Mesh topology enumeration + !! @param[in] coord_system Finite-element coordinate system enumeration + !! @param[in] scaled_radius Scaled planet radius + !! @param[in] rehabilitate Apply rehabilitation !! @param[in] ndf_w3 Number of degrees of freedom per cell for the operator space !! @param[in] basis_w3 Scalar basis functions evaluated at quadrature points !! @param[in] ndf_chi Number of degrees of freedom per cell for the coordinate field @@ -88,9 +96,12 @@ module sci_compute_mass_matrix_kernel_w3_mod ! R_SINGLE PRECISION ! ================== - subroutine compute_mass_matrix_w3_code_r_single( & + subroutine compute_mass_matrix_w3_code_r_single( & cell, nlayers, ncell_3d, mm, & chi1, chi2, chi3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & + rehabilitate, & ndf_w3, basis_w3, & ndf_chi, undf_chi, map_chi, & basis_chi, diff_basis_chi, & @@ -120,6 +131,12 @@ subroutine compute_mass_matrix_w3_code_r_single( & real(kind=r_single), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_single), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + logical(kind=l_def), intent(in) :: rehabilitate + !Internal variables integer(kind=i_def) :: df, df2, k, ik, ipanel integer(kind=i_def) :: qp1, qp2 @@ -184,6 +201,9 @@ end subroutine compute_mass_matrix_w3_code_r_single subroutine compute_mass_matrix_w3_code_mixed_precision( & cell, nlayers, ncell_3d, mm, & chi1, chi2, chi3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & + rehabilitate, & ndf_w3, basis_w3, & ndf_chi, undf_chi, map_chi, & basis_chi, diff_basis_chi, & @@ -213,6 +233,12 @@ subroutine compute_mass_matrix_w3_code_mixed_precision( & real(kind=r_double), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_double), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + logical, intent(in) :: rehabilitate + !Internal variables integer(kind=i_def) :: df, df2, k, ik, ipanel integer(kind=i_def) :: qp1, qp2 @@ -277,6 +303,9 @@ end subroutine compute_mass_matrix_w3_code_mixed_precision subroutine compute_mass_matrix_w3_code_r_double( & cell, nlayers, ncell_3d, mm, & chi1, chi2, chi3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & + rehabilitate, & ndf_w3, basis_w3, & ndf_chi, undf_chi, map_chi, & basis_chi, diff_basis_chi, & @@ -306,6 +335,12 @@ subroutine compute_mass_matrix_w3_code_r_double( & real(kind=r_double), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_double), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + logical, intent(in) :: rehabilitate + !Internal variables integer(kind=i_def) :: df, df2, k, ik, ipanel integer(kind=i_def) :: qp1, qp2 diff --git a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w_scalar_mod.F90 b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w_scalar_mod.F90 index 346ab1731..e020abafd 100644 --- a/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w_scalar_mod.F90 +++ b/components/science/source/kernel/fem/sci_compute_mass_matrix_kernel_w_scalar_mod.F90 @@ -15,20 +15,16 @@ module sci_compute_mass_matrix_kernel_w_scalar_mod GH_OPERATOR, GH_FIELD, & GH_LOGICAL, GH_SCALAR, & GH_READ, GH_WRITE, & - GH_REAL, ANY_SPACE_2, & - ANY_SPACE_9, & + GH_REAL, GH_INTEGER, & + ANY_SPACE_2, ANY_SPACE_9, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ - use constants_mod, only: i_def, r_single, r_double, l_def + use constants_mod, only: i_def, r_def, r_single, r_double, l_def use sci_coordinate_jacobian_mod, only: coordinate_jacobian use fs_continuity_mod, only: W0, Wtheta use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none private @@ -38,11 +34,15 @@ module sci_compute_mass_matrix_kernel_w_scalar_mod !--------------------------------------------------------------------------- type, public, extends(kernel_type) :: compute_mass_matrix_kernel_w_scalar_type private - type(arg_type) :: meta_args(4) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, ANY_SPACE_2, ANY_SPACE_2), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & - arg_type(GH_SCALAR, GH_LOGICAL, GH_READ) & + type(arg_type) :: meta_args(8) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, ANY_SPACE_2, ANY_SPACE_2), & ! mm + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_LOGICAL, GH_READ), & ! extended_mesh + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(2) = (/ & func_type(ANY_SPACE_2, GH_BASIS), & @@ -76,6 +76,10 @@ module sci_compute_mass_matrix_kernel_w_scalar_mod !! @param[in] chi3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels. !! @param[in] extended_mesh Compute on an extended mesh + !! @param[in] geometry Mesh geometry enumeration + !! @param[in] topology Mesh topology enumeration + !! @param[in] coord_system Finite-element coordinate system enumeration + !! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf_w_scalar The number of degrees of freedom per cell for w_scalar. !! @param[in] ndf_chi The number of degrees of freedom per cell for chi. !! @param[in] basis_w_scalar 4-dim array holding SCALAR basis functions evaluated at @@ -101,6 +105,8 @@ subroutine compute_mass_matrix_w_scalar_code_r32( & cell, nlayers, ncell_3d, mm, & chi1, chi2, chi3, panel_id, & extended_mesh, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w_scalar, basis_w_scalar, & ndf_chi, undf_chi, map_chi, & basis_chi, diff_basis_chi, & @@ -121,6 +127,11 @@ subroutine compute_mass_matrix_w_scalar_code_r32( & logical(kind=l_def), intent(in) :: extended_mesh + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + real(kind=r_single), dimension(ncell_3d,ndf_w_scalar,ndf_w_scalar), & intent(inout) :: mm @@ -201,6 +212,8 @@ subroutine compute_mass_matrix_w_scalar_code_r32r64( & cell, nlayers, ncell_3d, mm, & chi1, chi2, chi3, panel_id, & extended_mesh, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w_scalar, basis_w_scalar, & ndf_chi, undf_chi, map_chi, & basis_chi, diff_basis_chi, & @@ -221,6 +234,11 @@ subroutine compute_mass_matrix_w_scalar_code_r32r64( & logical(kind=l_def), intent(in) :: extended_mesh + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + real(kind=r_single), dimension(ncell_3d,ndf_w_scalar,ndf_w_scalar), & intent(inout) :: mm @@ -300,6 +318,8 @@ subroutine compute_mass_matrix_w_scalar_code_r64( & cell, nlayers, ncell_3d, mm, & chi1, chi2, chi3, panel_id, & extended_mesh, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w_scalar, basis_w_scalar, & ndf_chi, undf_chi, map_chi, & basis_chi, diff_basis_chi, & @@ -320,6 +340,11 @@ subroutine compute_mass_matrix_w_scalar_code_r64( & logical(kind=l_def), intent(in) :: extended_mesh + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + real(kind=r_double), dimension(ncell_3d,ndf_w_scalar,ndf_w_scalar), & intent(inout) :: mm diff --git a/components/science/source/kernel/fem/sci_gp_rhs_kernel_mod.F90 b/components/science/source/kernel/fem/sci_gp_rhs_kernel_mod.F90 index ca6d73d1e..9323af616 100644 --- a/components/science/source/kernel/fem/sci_gp_rhs_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_gp_rhs_kernel_mod.F90 @@ -7,14 +7,16 @@ !> @brief Kernel which projects a field into into a given space module sci_gp_rhs_kernel_mod + use kernel_mod, only : kernel_type use constants_mod, only : r_def, i_def -use argument_mod, only : arg_type, func_type, & - GH_FIELD, GH_REAL, GH_INC, & - GH_READ, ANY_SPACE_9, & - ANY_SPACE_1, ANY_SPACE_2, & - ANY_DISCONTINUOUS_SPACE_3, & - GH_BASIS, GH_DIFF_BASIS, & +use argument_mod, only : arg_type, func_type, & + GH_FIELD, GH_SCALAR, & + GH_REAL, GH_INTEGER, & + GH_INC, GH_READ, ANY_SPACE_9, & + ANY_SPACE_1, ANY_SPACE_2, & + ANY_DISCONTINUOUS_SPACE_3, & + GH_BASIS, GH_DIFF_BASIS, & CELL_COLUMN, GH_QUADRATURE_XYoZ implicit none @@ -27,11 +29,15 @@ module sci_gp_rhs_kernel_mod !> The type declaration for the kernel. Contains the metadata needed by the Psy layer type, public, extends(kernel_type) :: gp_rhs_kernel_type private - type(arg_type) :: meta_args(4) = (/ & - arg_type(GH_FIELD, GH_REAL, GH_INC, ANY_SPACE_1), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_SPACE_2), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(8) = (/ & + arg_type(GH_FIELD, GH_REAL, GH_INC, ANY_SPACE_1), & ! rhs + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_SPACE_2), & ! field + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(3) = (/ & func_type(ANY_SPACE_1, GH_BASIS), & @@ -63,6 +69,10 @@ module sci_gp_rhs_kernel_mod !! @param[in] chi_2 2nd coordinate field in Wchi !! @param[in] chi_3 3rd coordinate field in Wchi !! @param[in] panel_id Field giving the ID for mesh panels +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf Number of degrees of freedom per cell !! @param[in] undf Number of (local) unique degrees of freedom of the field rhs !! @param[in] map Dofmap for the cell at the base of the column @@ -87,6 +97,8 @@ module sci_gp_rhs_kernel_mod subroutine gp_rhs_code(nlayers, & rhs, field, & chi_1, chi_2, chi_3, panel_id, & + geometry, topology, & + coord_system, scaled_radius, & ndf, undf, map, basis, & ndf_f, undf_f, map_f, f_basis, & ndf_chi, undf_chi, map_chi, & @@ -96,10 +108,6 @@ subroutine gp_rhs_code(nlayers, & use sci_coordinate_jacobian_mod, only: coordinate_jacobian - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none ! Arguments @@ -124,6 +132,11 @@ subroutine gp_rhs_code(nlayers, & real(kind=r_def), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_def), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df2, k, qp1, qp2, ipanel real(kind=r_def), dimension(nqp_h,nqp_v) :: dj diff --git a/components/science/source/kernel/fem/sci_gp_vector_rhs_kernel_mod.F90 b/components/science/source/kernel/fem/sci_gp_vector_rhs_kernel_mod.F90 index 38a827c21..2b31d7d64 100644 --- a/components/science/source/kernel/fem/sci_gp_vector_rhs_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_gp_vector_rhs_kernel_mod.F90 @@ -8,8 +8,10 @@ module sci_gp_vector_rhs_kernel_mod use argument_mod, only : arg_type, func_type, & - GH_FIELD, GH_REAL, GH_INC, & - GH_READ, ANY_SPACE_1, & + GH_FIELD, GH_SCALAR, & + GH_REAL, GH_INTEGER, & + GH_INC, GH_READ, & + ANY_SPACE_1, & ANY_SPACE_2, ANY_SPACE_9, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & @@ -21,11 +23,7 @@ module sci_gp_vector_rhs_kernel_mod use coord_transform_mod, only : cart2sphere_vector use fs_continuity_mod, only : W0, W2 use kernel_mod, only : kernel_type - - use base_mesh_config_mod, only: geometry, topology, & - geometry_spherical - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius + use sci_mesh_enums_mod, only : geometry_spherical implicit none @@ -39,12 +37,16 @@ module sci_gp_vector_rhs_kernel_mod !> type, public, extends(kernel_type) :: gp_vector_rhs_kernel_type private - type(arg_type) :: meta_args(5) = (/ & - arg_type(GH_FIELD*3, GH_REAL, GH_INC, ANY_SPACE_1), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_SPACE_2), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & - arg_type(GH_FIELD, GH_REAL, GH_READ, W2) & + type(arg_type) :: meta_args(9) = (/ & + arg_type(GH_FIELD*3, GH_REAL, GH_INC, ANY_SPACE_1), & ! rhs1, rhs2, rhs3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_SPACE_2), & ! field + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_FIELD, GH_REAL, GH_READ, W2), & ! w2_field + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(3) = (/ & func_type(ANY_SPACE_1, GH_BASIS), & @@ -82,6 +84,10 @@ module sci_gp_vector_rhs_kernel_mod !! @param[in] chi_3 3rd coordinate field !! @param[in] panel_id Field giving the ID for mesh panels. !! @param[in] w2_field W2_field needed to get function space components +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf Number of degrees of freedom per cell !! @param[in] undf Number of degrees of freedom !! @param[in] map Dofmap for the cell at the base of the column @@ -111,6 +117,8 @@ subroutine gp_vector_rhs_code(nlayers, & rhs1, rhs2, rhs3, field, & chi_1, chi_2, chi_3, & panel_id, w2_field, & + geometry, topology, coord_system, & + scaled_radius, & ndf, undf, map, basis, & ndf_f, undf_f, map_f, f_basis, & ndf_chi, undf_chi, & @@ -149,6 +157,11 @@ subroutine gp_vector_rhs_code(nlayers, & real(kind=r_def), dimension(nqp_h), intent(in) :: wqp_h real(kind=r_def), dimension(nqp_v), intent(in) :: wqp_v + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df2, k, qp1, qp2 real(kind=r_def), dimension(nqp_h,nqp_v) :: dj diff --git a/components/science/source/kernel/fem/sci_mg_derham_mat_kernel_mod.F90 b/components/science/source/kernel/fem/sci_mg_derham_mat_kernel_mod.F90 index 04329bd5b..3143289f2 100644 --- a/components/science/source/kernel/fem/sci_mg_derham_mat_kernel_mod.F90 +++ b/components/science/source/kernel/fem/sci_mg_derham_mat_kernel_mod.F90 @@ -15,8 +15,10 @@ module sci_mg_derham_mat_kernel_mod use argument_mod, only: arg_type, func_type, & - GH_OPERATOR, GH_FIELD, & - GH_REAL, GH_READ, GH_WRITE, & + GH_OPERATOR, & + GH_FIELD, GH_SCALAR, & + GH_REAL, GH_INTEGER, & + GH_READ, GH_WRITE, & ANY_SPACE_1, ANY_SPACE_9, & ANY_DISCONTINUOUS_SPACE_3, & GH_BASIS, GH_DIFF_BASIS, & @@ -27,10 +29,6 @@ module sci_mg_derham_mat_kernel_mod use fs_continuity_mod, only: W2, W3, wtheta use kernel_mod, only: kernel_type - use base_mesh_config_mod, only: geometry, topology - use finite_element_config_mod, only: coord_system - use planet_config_mod, only: scaled_radius - implicit none private @@ -41,13 +39,17 @@ module sci_mg_derham_mat_kernel_mod type, public, extends(kernel_type) :: mg_derham_mat_kernel_type private - type(arg_type) :: meta_args(6) = (/ & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W2), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W3), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, Wtheta, Wtheta), & - arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2), & - arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & - arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3) & + type(arg_type) :: meta_args(10) = (/ & + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W2, W2), & ! mm2 + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W3), & ! mm3 + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, Wtheta, Wtheta), & ! mmt + arg_type(GH_OPERATOR, GH_REAL, GH_WRITE, W3, W2), & ! div + arg_type(GH_FIELD*3, GH_REAL, GH_READ, ANY_SPACE_9), & ! chi1, chi2, chi3 + arg_type(GH_FIELD, GH_REAL, GH_READ, ANY_DISCONTINUOUS_SPACE_3), & ! panel_id + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! geometry + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! topology + arg_type(GH_SCALAR, GH_INTEGER, GH_READ), & ! coord_system + arg_type(GH_SCALAR, GH_REAL, GH_READ) & ! scaled_radius /) type(func_type) :: meta_funcs(4) = (/ & func_type(W2, GH_BASIS, GH_DIFF_BASIS), & @@ -83,6 +85,10 @@ module sci_mg_derham_mat_kernel_mod !! @param[in] chi2 Physical coordinates in the 2nd dir !! @param[in] chi3 Physical coordinates in the 3rd dir !! @param[in] panel_id Field giving the ID for mesh panels +!! @param[in] geometry Mesh geometry enumeration +!! @param[in] topology Mesh topology enumeration +!! @param[in] coord_system Finite-element coordinate system enumeration +!! @param[in] scaled_radius Scaled planet radius !! @param[in] ndf_w2 Number of degrees of freedom per cell for W2 space !! @param[in] basis_w2 Basis functions evaluated at quadrature points for W2 space !! @param[in] diff_basis_w2 Differential of basis functions evaluated at @@ -112,6 +118,8 @@ subroutine mg_derham_mat_code(cell, nlayers, & ncell_3d6, div, & chi1, chi2, chi3, & panel_id, & + geometry, topology, & + coord_system, scaled_radius, & ndf_w2, basis_w2, diff_basis_w2, & ndf_w3, basis_w3, & ndf_wt, basis_wt, & @@ -149,6 +157,11 @@ subroutine mg_derham_mat_code(cell, nlayers, & real(kind=r_def), intent(in) :: wqp_h(nqp_h) real(kind=r_def), intent(in) :: wqp_v(nqp_v) + integer(kind=i_def), intent(in) :: geometry + integer(kind=i_def), intent(in) :: topology + integer(kind=i_def), intent(in) :: coord_system + real(kind=r_def), intent(in) :: scaled_radius + ! Internal variables integer(kind=i_def) :: df, df2, k, ik integer(kind=i_def) :: qp1, qp2 diff --git a/components/science/source/kernel/geometry/sci_chi_transform_mod.F90 b/components/science/source/kernel/geometry/sci_chi_transform_mod.F90 index 22e5bd04b..caad47efb 100644 --- a/components/science/source/kernel/geometry/sci_chi_transform_mod.F90 +++ b/components/science/source/kernel/geometry/sci_chi_transform_mod.F90 @@ -30,10 +30,11 @@ module sci_chi_transform_mod LOG_LEVEL_WARNING use matrix_invert_mod, only : matrix_invert_3x3 +use sci_mesh_enums_mod, only: geometry_spherical, & + geometry_planar, & + topology_fully_periodic + ! Configuration modules -use base_mesh_config_mod, only: geometry_spherical, & - geometry_planar, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_xyz, & coord_system_native diff --git a/components/science/source/kernel/geometry/sci_coordinate_jacobian_mod.F90 b/components/science/source/kernel/geometry/sci_coordinate_jacobian_mod.F90 index c0c0d4524..15ee48564 100644 --- a/components/science/source/kernel/geometry/sci_coordinate_jacobian_mod.F90 +++ b/components/science/source/kernel/geometry/sci_coordinate_jacobian_mod.F90 @@ -22,10 +22,10 @@ module sci_coordinate_jacobian_mod get_to_stretch, & get_to_rotate, & get_stretch_factor + use sci_mesh_enums_mod, only: geometry_planar, & + topology_fully_periodic ! Configuration modules - use base_mesh_config_mod, only: geometry_planar, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_xyz, & coord_system_native diff --git a/components/science/source/kernel/geometry/sci_height_continuous_kernel_mod.F90 b/components/science/source/kernel/geometry/sci_height_continuous_kernel_mod.F90 index 4964cb6aa..247a2e8a2 100644 --- a/components/science/source/kernel/geometry/sci_height_continuous_kernel_mod.F90 +++ b/components/science/source/kernel/geometry/sci_height_continuous_kernel_mod.F90 @@ -15,10 +15,10 @@ module sci_height_continuous_kernel_mod GH_READ, GH_INC, & ANY_SPACE_1, ANY_SPACE_9, & CELL_COLUMN, GH_BASIS, GH_EVALUATOR - use base_mesh_config_mod, only: geometry_spherical use constants_mod, only: r_def, i_def, l_def use finite_element_config_mod, only: coord_system_xyz use kernel_mod, only: kernel_type + use sci_mesh_enums_mod, only: geometry_spherical implicit none private diff --git a/components/science/source/kernel/geometry/sci_height_discontinuous_kernel_mod.F90 b/components/science/source/kernel/geometry/sci_height_discontinuous_kernel_mod.F90 index 101e8aa4b..f25dd4c6e 100644 --- a/components/science/source/kernel/geometry/sci_height_discontinuous_kernel_mod.F90 +++ b/components/science/source/kernel/geometry/sci_height_discontinuous_kernel_mod.F90 @@ -16,10 +16,10 @@ module sci_height_discontinuous_kernel_mod GH_READ, GH_WRITE, & ANY_DISCONTINUOUS_SPACE_1, ANY_SPACE_9, & CELL_COLUMN, GH_BASIS, GH_EVALUATOR - use base_mesh_config_mod, only: geometry_spherical use constants_mod, only: r_def, i_def, l_def use finite_element_config_mod, only: coord_system_xyz use kernel_mod, only: kernel_type + use sci_mesh_enums_mod, only: geometry_spherical implicit none private diff --git a/components/science/source/kernel/geometry/sci_native_jacobian_mod.F90 b/components/science/source/kernel/geometry/sci_native_jacobian_mod.F90 index 9e6b7446e..bec2de68b 100644 --- a/components/science/source/kernel/geometry/sci_native_jacobian_mod.F90 +++ b/components/science/source/kernel/geometry/sci_native_jacobian_mod.F90 @@ -24,11 +24,11 @@ module sci_native_jacobian_mod get_to_stretch, & get_to_rotate, & get_stretch_factor + use sci_mesh_enums_mod, only: geometry_planar, & + topology_fully_periodic use finite_element_config_mod, only: coord_system_xyz, & coord_system_native - use base_mesh_config_mod, only: geometry_planar, & - topology_fully_periodic implicit none diff --git a/components/science/source/kernel/inter_function_space/sci_compute_map_u_operators_kernel_mod.F90 b/components/science/source/kernel/inter_function_space/sci_compute_map_u_operators_kernel_mod.F90 index f05d600f8..dbaddf300 100644 --- a/components/science/source/kernel/inter_function_space/sci_compute_map_u_operators_kernel_mod.F90 +++ b/components/science/source/kernel/inter_function_space/sci_compute_map_u_operators_kernel_mod.F90 @@ -28,8 +28,7 @@ module sci_compute_map_u_operators_kernel_mod use fs_continuity_mod, only : W2, W3, Wtheta use kernel_mod, only : kernel_type use log_mod, only : log_event, LOG_LEVEL_ERROR, LOG_LEVEL_INFO - - use base_mesh_config_mod, only: geometry_spherical, geometry_planar + use sci_mesh_enums_mod, only : geometry_spherical, geometry_planar implicit none diff --git a/components/science/source/kernel/inter_function_space/sci_compute_sample_u_ops_kernel_mod.F90 b/components/science/source/kernel/inter_function_space/sci_compute_sample_u_ops_kernel_mod.F90 index eb8514b3c..9a1dc06f3 100644 --- a/components/science/source/kernel/inter_function_space/sci_compute_sample_u_ops_kernel_mod.F90 +++ b/components/science/source/kernel/inter_function_space/sci_compute_sample_u_ops_kernel_mod.F90 @@ -28,14 +28,14 @@ module sci_compute_sample_u_ops_kernel_mod use constants_mod, only : r_def, i_def use fs_continuity_mod, only : W2broken, W3, Wtheta use kernel_mod, only : kernel_type + use sci_chi_transform_mod, only : chi2llr use sci_coordinate_jacobian_mod, only : coordinate_jacobian, & coordinate_jacobian_inverse + use sci_mesh_enums_mod, only : geometry_spherical, geometry_planar use coord_transform_mod, only : sphere2cart_vector use reference_element_mod, only : W, S, N, E, T, B - use base_mesh_config_mod, only: geometry_spherical, geometry_planar - implicit none private diff --git a/components/science/source/kernel/inter_function_space/sci_convert_phys_to_hdiv_kernel_mod.F90 b/components/science/source/kernel/inter_function_space/sci_convert_phys_to_hdiv_kernel_mod.F90 index 2888e3577..9db530d0d 100644 --- a/components/science/source/kernel/inter_function_space/sci_convert_phys_to_hdiv_kernel_mod.F90 +++ b/components/science/source/kernel/inter_function_space/sci_convert_phys_to_hdiv_kernel_mod.F90 @@ -21,8 +21,7 @@ module sci_convert_phys_to_hdiv_kernel_mod use constants_mod, only : r_def, i_def use fs_continuity_mod, only : W2 use kernel_mod, only : kernel_type - - use base_mesh_config_mod, only: geometry_spherical + use sci_mesh_enums_mod, only : geometry_spherical implicit none diff --git a/components/science/source/kernel/inter_function_space/sci_project_ws_to_w1_operator_kernel_mod.F90 b/components/science/source/kernel/inter_function_space/sci_project_ws_to_w1_operator_kernel_mod.F90 index d00964b45..d8e569f88 100644 --- a/components/science/source/kernel/inter_function_space/sci_project_ws_to_w1_operator_kernel_mod.F90 +++ b/components/science/source/kernel/inter_function_space/sci_project_ws_to_w1_operator_kernel_mod.F90 @@ -26,8 +26,7 @@ module sci_project_ws_to_w1_operator_kernel_mod use constants_mod, only : r_def, i_def use fs_continuity_mod, only : W1 use log_mod, only : log_event, LOG_LEVEL_ERROR - -use base_mesh_config_mod, only: geometry_spherical, geometry_planar +use sci_mesh_enums_mod, only : geometry_spherical, geometry_planar implicit none diff --git a/components/science/source/solver/sci_hierarchical_linear_operator_mod.f90 b/components/science/source/solver/sci_hierarchical_linear_operator_mod.f90 index 69cb77928..39dd957bc 100644 --- a/components/science/source/solver/sci_hierarchical_linear_operator_mod.f90 +++ b/components/science/source/solver/sci_hierarchical_linear_operator_mod.f90 @@ -12,6 +12,7 @@ module sci_hierarchical_linear_operator_mod use sci_linear_operator_mod, only : abstract_linear_operator_type + use config_mod, only: config_type implicit none @@ -28,12 +29,15 @@ module sci_hierarchical_linear_operator_mod !> Abstract interface defined for the coarsening linear operator to the !> next level of the multigrid hierarchy. !> - !>@param[in] self a hierarchical linear operator - !>@param[in] fs_coarse coarse level function space - !>@param[inout] coarse_operator coarsened version on the next multigrid level - subroutine coarsen_interface(self,coarse_operator) + !>@param[in] config Application configuration object + !>@param[inout] coarse_operator Coarsened version on the next multigrid level + subroutine coarsen_interface(self, config, coarse_operator) import :: abstract_hierarchical_linear_operator_type + import :: config_type + class(abstract_hierarchical_linear_operator_type), intent(in) :: self + + type(config_type), intent(in) :: config class(abstract_hierarchical_linear_operator_type), allocatable, intent(inout) :: coarse_operator end subroutine coarsen_interface end interface diff --git a/components/science/source/solver/sci_hierarchical_preconditioner_mod.f90 b/components/science/source/solver/sci_hierarchical_preconditioner_mod.f90 index e63d8f42f..040be6cf7 100644 --- a/components/science/source/solver/sci_hierarchical_preconditioner_mod.f90 +++ b/components/science/source/solver/sci_hierarchical_preconditioner_mod.f90 @@ -14,6 +14,7 @@ module sci_hierarchical_preconditioner_mod use sci_preconditioner_mod, only: abstract_preconditioner_type + use config_mod, only: config_type implicit none private @@ -29,11 +30,15 @@ module sci_hierarchical_preconditioner_mod !> Abstract interface defined for the coarsening preconditioner to the !> next level of the multigrid hierarchy. !> - !>@param[inout] self a hierarchical preconditioner - !>@param[inout] other coarsened version on the next multigrid level - subroutine coarsen_interface(self, other) + !>@param[in] config Application configuration object + !>@param[inout] other Coarsened version on the next multigrid level + subroutine coarsen_interface(self, config, other) import :: abstract_hierarchical_preconditioner_type + import :: config_type + class(abstract_hierarchical_preconditioner_type), intent(inout) :: self + + type(config_type), intent(in) :: config class(abstract_hierarchical_preconditioner_type), allocatable, intent(inout) :: other end subroutine coarsen_interface end interface diff --git a/components/science/source/solver/sci_iterative_solver_mod.f90 b/components/science/source/solver/sci_iterative_solver_mod.f90 index f85ef074f..b955242f2 100644 --- a/components/science/source/solver/sci_iterative_solver_mod.f90 +++ b/components/science/source/solver/sci_iterative_solver_mod.f90 @@ -14,6 +14,7 @@ module sci_iterative_solver_mod + use config_mod, only : config_type use constants_mod, only : r_def, i_def, EPS, l_def use vector_mod, only : abstract_vector_type use sci_linear_operator_mod, & @@ -63,15 +64,22 @@ module sci_iterative_solver_mod !> @details Apply the iterative solver for a given right hand side !> \f$b\f$. !> - !> @param[inout] x Resulting solution \f$x\f$ - !> @param[in] b Right hand side vector \f$b\f$ + !> @param[in] config Application configuration object + !> @param[inout] x Resulting solution \f$x\f$ + !> @param[in] b Right hand side vector \f$b\f$ !> - subroutine apply_interface(self, x, b) + subroutine apply_interface(self, config, x, b) + import :: abstract_vector_type import :: abstract_iterative_solver_type + import :: config_type + class(abstract_iterative_solver_type), intent(inout) :: self + + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b + end subroutine apply_interface end interface @@ -113,8 +121,9 @@ module function cg_constructor( lin_op, prec, r_tol, a_tol, max_iter, & end interface interface - module subroutine cg_solve(self, x, b) + module subroutine cg_solve(self, config, x, b) class(conjugate_gradient_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine @@ -149,8 +158,9 @@ module function bicgstab_constructor( lin_op, prec, r_tol, a_tol, max_iter, & end interface interface - module subroutine bicgstab_solve(self, x, b) + module subroutine bicgstab_solve(self, config, x, b) class(bicgstab_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine @@ -188,8 +198,9 @@ end function gmres_constructor end interface interface - module subroutine gmres_solve(self, x, b) + module subroutine gmres_solve(self, config, x, b) class(gmres_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine gmres_solve @@ -226,8 +237,9 @@ end function fgmres_constructor end interface interface - module subroutine fgmres_solve(self, x, b) - class(fgmres_type), intent(inout) :: self + module subroutine fgmres_solve(self, config, x, b) + class(fgmres_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine fgmres_solve @@ -264,8 +276,9 @@ end function gcr_constructor end interface interface - module subroutine gcr_solve(self, x, b) + module subroutine gcr_solve(self, config, x, b) class(gcr_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine gcr_solve @@ -304,8 +317,9 @@ end function block_gcr_constructor end interface interface - module subroutine block_gcr_solve(self, x, b) + module subroutine block_gcr_solve(self, config, x, b) class(block_gcr_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine block_gcr_solve @@ -339,8 +353,9 @@ module function precondition_only_constructor( lin_op, prec, & end interface interface - module subroutine precondition_only_solve(self, x, b) + module subroutine precondition_only_solve(self, config, x, b) class(precondition_only_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine @@ -380,10 +395,11 @@ module function jacobi_constructor( lin_op, prec, r_tol, a_tol, max_iter, end interface interface - module subroutine jacobi_solve(self, x, b) - class(jacobi_type), intent(inout) :: self - class(abstract_vector_type), intent(inout) :: x - class(abstract_vector_type), intent(inout) :: b + module subroutine jacobi_solve(self, config, x, b) + class(jacobi_type), intent(inout) :: self + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(inout) :: x + class(abstract_vector_type), intent(inout) :: b end subroutine end interface @@ -421,8 +437,9 @@ module function chebyshev_constructor( lin_op, prec, r_tol, a_tol, max_iter, & end interface interface - module subroutine chebyshev_solve(self, x, b) + module subroutine chebyshev_solve(self, config, x, b) class(chebyshev_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b end subroutine @@ -490,12 +507,16 @@ module function cg_constructor(lin_op, prec, r_tol, a_tol, max_iter, & !> !> Over-rides the abstract interface to do the actual solve. !> - !> @param[inout] b "RHS" or boundary conditions. - !> @param[inout] x Solution. + !> @param[in] config Application configuration object + !> @param[inout] b "RHS" or boundary conditions. + !> @param[inout] x Solution. !> - module subroutine cg_solve(self, x, b) + module subroutine cg_solve(self, config, x, b) + implicit none + class(conjugate_gradient_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b @@ -519,7 +540,7 @@ module subroutine cg_solve(self, x, b) converged=.false. !set up the algorithm - call self%lin_op%apply(x,r) ! r = A.x + call self%lin_op%apply(config, x, r) ! r = A.x call r%scale(-1.0_r_def) ! r = -A.x call r%axpy(1.0_r_def, b) ! r = b - A.x r_nrm_0 = r%norm() ! r_0 = ||r||_2 @@ -539,8 +560,8 @@ module subroutine cg_solve(self, x, b) end if call z%set_scalar(0.0_r_def) - call self%prec%apply(r,z) ! z = P^{-1}.r - rz = r%dot(z) ! rz = + call self%prec%apply(config, r, z) ! z = P^{-1}.r + rz = r%dot(z) ! rz = r_nrm_old = r_nrm_0 call p%copy(z) @@ -550,10 +571,10 @@ module subroutine cg_solve(self, x, b) end if ! iterate until maximal number of iterations is reached do iter=1, self%max_iter - call self%lin_op%apply(p,z) ! z = A.p - alpha = rz / p%dot(z) ! alpha = / - call x%axpy(alpha,p) ! x -> x + alpha*p - call r%axpy(-alpha,z) ! r -> r - alpha*A.p + call self%lin_op%apply(config, p, z) ! z = A.p + alpha = rz / p%dot(z) ! alpha = / + call x%axpy(alpha,p) ! x -> x + alpha*p + call r%axpy(-alpha,z) ! r -> r - alpha*A.p if ( self%monitor_convergence ) then r_nrm = r%norm() ! r = ||r||_2 @@ -567,10 +588,10 @@ module subroutine cg_solve(self, x, b) exit end if end if - call self%prec%apply(r,z) ! z = P^{-1}.r - rz_new = r%dot(z) ! rz_new = - beta = rz_new/rz ! beta = / - call p%aypx(beta,z) ! p -> z + beta*p + call self%prec%apply(config, r, z) ! z = P^{-1}.r + rz_new = r%dot(z) ! rz_new = + beta = rz_new/rz ! beta = / + call p%aypx(beta,z) ! p -> z + beta*p rz = rz_new r_nrm_old = r_nrm end do @@ -640,15 +661,19 @@ module function bicgstab_constructor( lin_op, prec, r_tol, a_tol, max_iter, & end function bicgstab_constructor !> bicgstab solve. Over-rides the abstract interface to do the actual solve. - !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type - !> This the "RHS" or boundary conditions, - !> @param[inout] x an abstract vector which is the solution + !> @param[in] config Application configuration object + !> @param[inout] b An abstract vector which will be an actual vector of + !> unkown extended type. This the "RHS" or boundary conditions, + !> @param[inout] x An abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine bicgstab_solve(self, x, b) + module subroutine bicgstab_solve(self, config, x, b) + implicit none - class(bicgstab_type), intent(inout) :: self - class(abstract_vector_type), intent(inout) :: x - class(abstract_vector_type), intent(inout) :: b + + class(bicgstab_type), intent(inout) :: self + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(inout) :: x + class(abstract_vector_type), intent(inout) :: b ! tempory vectors class(abstract_vector_type), allocatable :: r @@ -674,7 +699,7 @@ module subroutine bicgstab_solve(self, x, b) ! v = Ax call x%duplicate(v) call v%set_scalar(0.0_r_def) - call self%lin_op%apply(x,v) + call self%lin_op%apply(config, x, v) ! r = b - Ax call r%axpy(-1.0_r_def,v) ! store initial residual @@ -721,20 +746,20 @@ module subroutine bicgstab_solve(self, x, b) call t%copy(r) call t%axpy(-beta*omega, v) ! stage 2 post-condition - call self%prec%apply(t,y) + call self%prec%apply(config, t, y) ! now add on beta P call p%aypx(beta, y) ! apply the matrix - call self%lin_op%apply(p,v) + call self%lin_op%apply(config, p, v) alpha = rho/r0%dot(v) ! s = r - alpha * v call s%copy(r) call s%axpy(-alpha,v) ! apply the preconditioner - call self%prec%apply(s,z) + call self%prec%apply(config, s, z) ! apply the operator - call self%lin_op%apply(z,t) + call self%lin_op%apply(config, z, t) ! final scalars tt = t%dot(t) @@ -826,15 +851,19 @@ end function gmres_constructor !> gmres_solve. Over-rides the abstract interface to do the actual solve. !> @detail The solver implements left-preconditioning, i.e. is solving M{-1}.A.x = M{-1}.b + !> @param[in] config Application configuration object !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type !> This the "RHS" or boundary conditions, !> @param[inout] x an abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine gmres_solve(self, x, b) + module subroutine gmres_solve(self, config, x, b) + implicit none - class(gmres_type), intent(inout) :: self - class(abstract_vector_type), intent(inout) :: x - class(abstract_vector_type), intent(inout) :: b + + class(gmres_type), intent(inout) :: self + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(inout) :: x + class(abstract_vector_type), intent(inout) :: b ! temporary vectors class(abstract_vector_type), allocatable :: s @@ -859,7 +888,7 @@ module subroutine gmres_solve(self, x, b) call b%duplicate(res) ! compute res = b -Ax ... in stages - call self%lin_op%apply(x,Ax) + call self%lin_op%apply(config, x, Ax) call res%copy(b) call res%axpy(-1.0_r_def,Ax) @@ -899,7 +928,7 @@ module subroutine gmres_solve(self, x, b) ! initialisation complete, lets go to work. do iter = 1, self%max_iter call s%set_scalar(0.0_r_def) - call self%prec%apply(res,s) + call self%prec%apply(config, res, s) beta = s%norm() call v(1)%vt%copy(s) call v(1)%vt%scale(1.0_r_def/beta) @@ -911,9 +940,9 @@ module subroutine gmres_solve(self, x, b) call w%copy(v(iv)%vt) ! apply the operator - call self%lin_op%apply( w, s ) + call self%lin_op%apply(config, w, s) ! apply the preconditioner - call self%prec%apply( s, w ) + call self%prec%apply(config, s, w) ! compute the h values do ivj = 1, iv @@ -963,7 +992,7 @@ module subroutine gmres_solve(self, x, b) end do call Ax%set_scalar(0.0_r_def) - call self%lin_op%apply(x, Ax) + call self%lin_op%apply(config, x, Ax) call res%copy(Ax) call res%aypx(-1.0_r_def, b) @@ -1052,15 +1081,17 @@ end function fgmres_constructor !> @detail The solver implements flexible right-preconditioning, i.e. is solving A.M{-1}.M.x = b !> as 1) A.M{-1}.y = b !> 2) M{-1}y = x + !> @param[in] config Application configuration object !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type !> This the "RHS" or boundary conditions, !> @param[inout] x an abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine fgmres_solve(self, x, b) + module subroutine fgmres_solve(self, config, x, b) implicit none - class(fgmres_type), intent(inout) :: self - class(abstract_vector_type), intent(inout) :: x - class(abstract_vector_type), intent(inout) :: b + class(fgmres_type), intent(inout) :: self + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(inout) :: x + class(abstract_vector_type), intent(inout) :: b ! temporary vectors class(abstract_vector_type), allocatable :: s @@ -1138,9 +1169,9 @@ module subroutine fgmres_solve(self, x, b) do iter = 1, self%max_iter do iv = 1, self%gcrk - call self%prec%apply(v(iv)%vt, Pv(iv)%vt) + call self%prec%apply(config, v(iv)%vt, Pv(iv)%vt) ! apply the operator - call self%lin_op%apply( Pv(iv)%vt, s ) + call self%lin_op%apply(config, Pv(iv)%vt, s) call w%copy(s) ! compute the h values do ivj = 1, iv @@ -1189,7 +1220,7 @@ module subroutine fgmres_solve(self, x, b) end do ! check for convergence - call self%lin_op%apply(dx, Ax) + call self%lin_op%apply(config, dx, Ax) call res%copy(Ax) call res%aypx(-1.0_r_def, b) @@ -1286,15 +1317,19 @@ end function gcr_constructor !> @detail The solver implements right-preconditioning, i.e. is solving A.M{-1}.M.x = b !> as 1) A.M{-1}.y = b !> 2) M{-1}y = x + !> @param[in] config Application configuration object !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type !> This the "RHS" or boundary conditions, !> @param[inout] x an abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine gcr_solve(self, x, b) + module subroutine gcr_solve(self, config, x, b) + implicit none - class(gcr_type), intent(inout) :: self - class(abstract_vector_type), intent(inout) :: x - class(abstract_vector_type), intent(inout) :: b + + class(gcr_type), intent(inout) :: self + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(inout) :: x + class(abstract_vector_type), intent(inout) :: b ! temporary vectors class(abstract_vector_type), allocatable :: dx @@ -1317,8 +1352,8 @@ module subroutine gcr_solve(self, x, b) ! initial guess call dx%set_scalar(0.0_r_def) - call self%prec%apply( b, dx ) - call self%lin_op%apply( dx, Ax ) + call self%prec%apply(config, b, dx) + call self%lin_op%apply(config, dx, Ax) !res = b - Ax call res%copy(b) call res%axpy(-1.0_r_def, Ax) @@ -1356,9 +1391,9 @@ module subroutine gcr_solve(self, x, b) do iv = 1, self%gcrk ! apply the preconditioner - call self%prec%apply( res, Pv(iv)%vt ) + call self%prec%apply(config, res, Pv(iv)%vt) ! apply the operator - call self%lin_op%apply( Pv(iv)%vt, v(iv)%vt ) + call self%lin_op%apply(config, Pv(iv)%vt, v(iv)%vt) do ivj = 1, iv-1 alpha = v(iv)%vt%dot(v(ivj)%vt) @@ -1460,15 +1495,17 @@ module function block_gcr_constructor( lin_op, prec, gcrk, r_tol, a_tol, max_ite !> @detail The solver implements right-preconditioning, i.e. is solving A.M{-1}.M.x = b !> as 1) A.M{-1}.y = b !> 2) M{-1}y = x + !> @param[in] config Application configuration object !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type !> This the "RHS" or boundary conditions, !> @param[inout] x an abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine block_gcr_solve(self, x, b) + module subroutine block_gcr_solve(self, config, x, b) implicit none - class(block_gcr_type), intent(inout) :: self - class(abstract_vector_type), intent(inout) :: x - class(abstract_vector_type), intent(inout) :: b + class(block_gcr_type), intent(inout) :: self + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(inout) :: x + class(abstract_vector_type), intent(inout) :: b ! temporary vectors class(abstract_vector_type), allocatable :: Ax @@ -1501,8 +1538,8 @@ module subroutine block_gcr_solve(self, x, b) ! initial guess call x%set_scalar(0.0_r_def) - call self%prec%apply( b, x ) - call self%lin_op%apply( x, Ax ) + call self%prec%apply(config, b, x) + call self%lin_op%apply(config, x, Ax) call res%copy(b) call res%axpy(-1.0_r_def, Ax) @@ -1536,9 +1573,9 @@ module subroutine block_gcr_solve(self, x, b) do iv = 1, self%gcrk ! apply the preconditioner - call self%prec%apply( res, Pv(iv)%vt ) + call self%prec%apply(config, res, Pv(iv)%vt) ! apply the operator - call self%lin_op%apply( Pv(iv)%vt, v(iv)%vt ) + call self%lin_op%apply(config, Pv(iv)%vt, v(iv)%vt) do ivj = 1, iv-1 alpha = v(iv)%vt%dot(v(ivj)%vt) @@ -1661,13 +1698,17 @@ module function precondition_only_constructor(lin_op, prec, monitor_convergence) end function !> Precondition only solve. Over-rides the abstract interface to do the actual solve. + !> @param[in] config Application configuration object !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type !> This the "RHS" or boundary conditions, !> @param[inout] x an abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine precondition_only_solve(self, x, b) + module subroutine precondition_only_solve(self, config, x, b) + implicit none + class(precondition_only_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b @@ -1677,7 +1718,7 @@ module subroutine precondition_only_solve(self, x, b) call log_event("Precondition only starting", LOG_LEVEL_DEBUG) - call self%prec%apply(b,x) ! x = P^{-1}.b + call self%prec%apply(config, b, x) ! x = P^{-1}.b if( self%monitor_convergence ) then ! Compute initial and final error @@ -1685,7 +1726,7 @@ module subroutine precondition_only_solve(self, x, b) call res%copy(b) e0 = res%norm() call x%duplicate(Ax) - call self%lin_op%apply(x, Ax) + call self%lin_op%apply(config, x, Ax) call res%axpy(-1.0_r_def, Ax) e = res%norm() write(log_scratch_space,'(A,3E15.8)') & @@ -1743,13 +1784,17 @@ module function jacobi_constructor(lin_op, prec, r_tol, a_tol, max_iter, & end function !> Jacobi solve. Over-rides the abstract interface to do the actual solve. + !> @param[in] config Application configuration object !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type !> This the "RHS" or boundary conditions, !> @param[inout] x an abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine jacobi_solve(self, x, b) + module subroutine jacobi_solve(self, config, x, b) + implicit none + class(jacobi_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b @@ -1774,7 +1819,7 @@ module subroutine jacobi_solve(self, x, b) ! Solve Ax = b do iter=1, self%max_iter - call self%lin_op%apply(x, r) ! r = Ax + call self%lin_op%apply(config, x, r) ! r = Ax call r%axpy(-1.0_r_def, b) ! r = r - b = Ax - b ! Check for convergence @@ -1792,7 +1837,7 @@ module subroutine jacobi_solve(self, x, b) end if end if - call self%prec%apply(r, z) ! z = P^{-1}.r = P^{-1}.(Ax - b) + call self%prec%apply(config, r, z) ! z = P^{-1}.r = P^{-1}.(Ax - b) call x%axpy(const, z) ! x = x + const*z = x + const*P^{-1}(Ax - b) end do @@ -1867,13 +1912,17 @@ module function chebyshev_constructor(lin_op, prec, r_tol, a_tol, max_iter, end function !> chebyshev solve. Over-rides the abstract interface to do the actual solve. + !> @param[in] config Application configuration object !> @param[inout] b an abstract vector which will be an actual vector of unkown extended type !> This the "RHS" or boundary conditions, !> @param[inout] x an abstract vector which is the solution !> @param[self] The solver which has pointers to the lin_op and preconditioner - module subroutine chebyshev_solve(self, x, b) + module subroutine chebyshev_solve(self, config, x, b) + implicit none + class(chebyshev_type), intent(inout) :: self + type(config_type), intent(in) :: config class(abstract_vector_type), intent(inout) :: x class(abstract_vector_type), intent(inout) :: b @@ -1908,11 +1957,11 @@ module subroutine chebyshev_solve(self, x, b) do iter = 1, self%max_iter ! r = b-M*xo - call self%lin_op%apply(xo,z) + call self%lin_op%apply(config, xo, z) call r%axpby(1.0_r_def, b, -1.0_r_def, z) ! z = D^{-1}.r - call self%prec%apply(r,z) + call self%prec%apply(config, r, z) ! x = w*(a/b*z+xo) + (1-w)*xp w = 1.0_r_def/(1.0_r_def - w/(4.0_r_def*a2**2)) @@ -1931,7 +1980,7 @@ module subroutine chebyshev_solve(self, x, b) end do ! residiual = norm(b - M*x) if ( self%monitor_convergence ) then - call self%lin_op%apply(x,z) ! z = M.x + call self%lin_op%apply(config, x, z) ! z = M.x call z%axpy(-1.0_r_def, b) ! z = M.x-b final_norm = z%norm() write(log_scratch_space, & diff --git a/components/science/source/solver/sci_linear_operator_mod.f90 b/components/science/source/solver/sci_linear_operator_mod.f90 index 87ef31dab..78c8bef41 100644 --- a/components/science/source/solver/sci_linear_operator_mod.f90 +++ b/components/science/source/solver/sci_linear_operator_mod.f90 @@ -9,8 +9,11 @@ !! defines an interface for the linear operator application y = A.x module sci_linear_operator_mod - use vector_mod, only : abstract_vector_type - use function_space_mod, only : function_space_type + + use config_mod, only: config_type + use function_space_mod, only: function_space_type + use vector_mod, only: abstract_vector_type + implicit none private @@ -22,16 +25,21 @@ module sci_linear_operator_mod end type abstract_linear_operator_type abstract interface - !> Abstract interface defined for the apply method. - !> param[in] self a linear operator - !> param[in] x a vector the linear operator is applied to - !> param[inout] y a vector, the result. - subroutine apply_interface(self, x, y) + !> @brief Abstract interface defined for the apply method. + !> @param[in] config Application configuration object + !> @param[in] x a vector the linear operator is applied to + !> @param[inout] y a vector, the result. + subroutine apply_interface(self, config, x, y) import :: abstract_linear_operator_type import :: abstract_vector_type + import :: config_type + class(abstract_linear_operator_type), intent(inout) :: self - class(abstract_vector_type), intent(in) :: x - class(abstract_vector_type), intent(inout) :: y + + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(in) :: x + class(abstract_vector_type), intent(inout) :: y + end subroutine apply_interface end interface diff --git a/components/science/source/solver/sci_preconditioner_mod.f90 b/components/science/source/solver/sci_preconditioner_mod.f90 index 2b88c9a6a..a602a2961 100644 --- a/components/science/source/solver/sci_preconditioner_mod.f90 +++ b/components/science/source/solver/sci_preconditioner_mod.f90 @@ -10,7 +10,9 @@ !> defines an interface for the preconditioner application y = P^{-1}.x module sci_preconditioner_mod - use vector_mod, only : abstract_vector_type + + use config_mod, only: config_type + use vector_mod, only: abstract_vector_type implicit none private @@ -25,15 +27,20 @@ module sci_preconditioner_mod abstract interface !> abstract interface defined for the apply procedure of a preconditioner !! y = P.x - !> @param[in] self a preconditioner - !> @param[in] x a vector that the preconditioner is applied to. - !> @param[inout] y a vector, the result. - subroutine apply_interface(self, x, y) + !> @param[in] config Application configuration object + !> @param[in] x A vector that the preconditioner is applied to. + !> @param[inout] y A vector, the result. + subroutine apply_interface(self, config, x, y) + import :: abstract_vector_type import :: abstract_preconditioner_type + import :: config_type + class(abstract_preconditioner_type), intent(inout) :: self - class(abstract_vector_type), intent(in) :: x - class(abstract_vector_type), intent(inout) :: y + + type(config_type), intent(in) :: config + class(abstract_vector_type), intent(in) :: x + class(abstract_vector_type), intent(inout) :: y end subroutine apply_interface end interface diff --git a/components/science/source/utilities/sci_mesh_enums_mod.f90 b/components/science/source/utilities/sci_mesh_enums_mod.f90 new file mode 100644 index 000000000..a26d37ab6 --- /dev/null +++ b/components/science/source/utilities/sci_mesh_enums_mod.f90 @@ -0,0 +1,77 @@ +!----------------------------------------------------------------------------- +! (c) Crown copyright Met Office. All rights reserved. +! The file LICENCE, distributed with this code, contains details of the terms +! under which the code may be used +!----------------------------------------------------------------------------- + +!> @brief Module to access mesh enumerations. +module sci_mesh_enums_mod + + use constants_mod, only: i_def, imdi + use log_mod, only: log_event, log_level_error + use mesh_mod, only: mesh_type + + use base_mesh_config_mod, only: & + sci_geometry_spherical => geometry_spherical, & + sci_geometry_planar => geometry_planar, & + sci_topology_fully_periodic => topology_fully_periodic, & + sci_topology_non_periodic => topology_non_periodic + + implicit none + + private + + public :: geometry_spherical, geometry_planar + public :: topology_fully_periodic, topology_non_periodic + + public :: get_mesh_enums + + ! These will get switched to something hardcoded for the science components + ! at a later date to break dependence on base_mesh_config_mod + integer(i_def), parameter :: geometry_spherical = sci_geometry_spherical ! 157 + integer(i_def), parameter :: geometry_planar = sci_geometry_planar ! 358 + integer(i_def), parameter :: topology_fully_periodic = sci_topology_fully_periodic ! 492 + integer(i_def), parameter :: topology_non_periodic = sci_topology_non_periodic ! 157 + +contains + +!--------------------------------------------------------------------------- +!> @brief Returns mesh enumerations in line with science component values +!> @param[in] mesh Mesh object to query +!> @param[out] geometry [optional] Science component geometry enumeration +!> @param[out] topology [optional] Science component topology enumeration +!> +subroutine get_mesh_enums(mesh, geometry, topology) + + implicit none + + type(mesh_type), intent(in) :: mesh + + integer(i_def), optional, intent(out) :: geometry + integer(i_def), optional, intent(out) :: topology + + if (present(geometry)) then + geometry = imdi + if (mesh%is_geometry_spherical()) then + geometry = geometry_spherical + else if (mesh%is_geometry_planar()) then + geometry = geometry_planar + else + call log_event('Unsupported mesh geometry', log_level_error) + end if + end if + + if (present(topology)) then + topology = imdi + if (mesh%is_topology_periodic()) then + topology = topology_fully_periodic + else if (mesh%is_topology_non_periodic()) then + topology = topology_non_periodic + else + call log_event('Unsupported mesh topology', log_level_error) + end if + end if + +end subroutine get_mesh_enums + +end module sci_mesh_enums_mod diff --git a/components/science/source/utilities/sci_query_mod.f90 b/components/science/source/utilities/sci_query_mod.f90 index 7613aabcc..1db2491f2 100644 --- a/components/science/source/utilities/sci_query_mod.f90 +++ b/components/science/source/utilities/sci_query_mod.f90 @@ -16,9 +16,9 @@ module sci_query_mod implicit none private - public :: valid_for_global_model, & - check_lbc, & - is_lbc + public :: valid_for_global_model + public :: check_lbc + public :: is_lbc interface is_lbc module procedure is_global_lbc @@ -109,7 +109,6 @@ function check_lbc( mesh_name ) result ( answer ) character(4), parameter :: lbc_tag = '-lbc' integer(i_def) :: mesh_name_length -! character(4) :: tmp_str ! Check to see if the loaded mesh is an LBC mesh ! where the lbc_suffix would be appended. diff --git a/components/science/unit-test/kernel/fem/compute_broken_div_operator_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/compute_broken_div_operator_kernel_mod_test.pf index ea0f476ef..e22e9ce2d 100644 --- a/components/science/unit-test/kernel/fem/compute_broken_div_operator_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_broken_div_operator_kernel_mod_test.pf @@ -7,7 +7,9 @@ !> module compute_broken_div_operator_kernel_mod_test - use constants_mod, only : i_def, r_def, imdi + use constants_mod, only : l_def, i_def, r_def, imdi, rmdi + use finite_element_config_mod, only: coord_system_xyz + use funit use get_unit_test_m3x3_dofmap_mod, & only : get_w0_m3x3_dofmap, & get_w3_m3x3_dofmap @@ -25,76 +27,47 @@ module compute_broken_div_operator_kernel_mod_test get_gaussian_q3x3x3_quadrature_weights_z use get_unit_test_3x3x3_chi_mod, & only : get_w0_3x3x3_field - use funit implicit none private - public :: test_all - - @TestCase - type, extends(TestCase), public :: compute_broken_div_operator_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_broken_div_operator_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_broken_div_operator_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=1_i_def, & - element_order_v=1_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_broken_div_operator_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_all( this ) + subroutine test_all() use sci_compute_broken_div_operator_kernel_mod, & only : compute_broken_div_operator_code implicit none - class(compute_broken_div_operator_test_type), intent(inout) :: this - real(r_def), parameter :: tol = 1.0e-6_r_def real(r_def), parameter :: dx = 6000.0_r_def real(r_def), parameter :: dy = 1000.0_r_def @@ -102,7 +75,12 @@ contains integer(i_def), parameter :: nlayers = 3 - integer(i_def) :: i, j, k, icell + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + logical(l_def), parameter :: rehabilitate = .true. + integer(i_def) :: ndf_w0, undf_w0, ndf_w2b integer(i_def) :: ndf_w3, ndf_pid, undf_pid integer(i_def) :: nqp_h, nqp_v @@ -165,6 +143,11 @@ contains chi2, & chi3, & panel_id, & + geometry, & + topology, & + coord_system, & + scaled_radius, & + rehabilitate, & ndf_w3, & basis_w3, & ndf_w2b, & diff --git a/components/science/unit-test/kernel/fem/compute_curl_operator_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/compute_curl_operator_kernel_mod_test.pf index 12ae5c4c5..06b536476 100644 --- a/components/science/unit-test/kernel/fem/compute_curl_operator_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_curl_operator_kernel_mod_test.pf @@ -8,7 +8,7 @@ !> module compute_curl_operator_kernel_mod_test - use constants_mod, only : i_def, r_def, imdi + use constants_mod, only : i_def, r_def, imdi, rmdi use get_unit_test_m3x3_dofmap_mod, only : get_w0_m3x3_dofmap, & get_w3_m3x3_dofmap use get_unit_test_m3x3_q3x3x3_sizes_mod, only : get_w0_m3x3_q3x3x3_size, & @@ -25,87 +25,62 @@ module compute_curl_operator_kernel_mod_test use get_unit_test_3x3x3_chi_mod, only : get_w0_3x3x3_field use funit + use finite_element_config_mod, only: coord_system_xyz + implicit none private - public :: test_all - @TestCase - type, extends(TestCase), public :: compute_curl_operator_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_curl_operator_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_curl_operator_test_type), intent(inout) :: this - - integer(i_def) :: nlayers, dim, diff_dim, dim_w2, diff_dim_w1 - integer(i_def) :: ndf_w0, ndf_w2, ndf_w1 - integer(i_def) :: nqp_h, nqp_v - integer(i_def) :: i - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_curl_operator_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_all( this ) + subroutine test_all() use sci_compute_curl_operator_kernel_mod, & only : compute_curl_operator_code implicit none - class(compute_curl_operator_test_type), intent(inout) :: this - real(kind=r_def), parameter :: tol = 1.0e-6_r_def - real(kind=r_def), parameter :: dx = 6000.0_r_def - real(r_def), parameter :: dy = 1000.0_r_def - real(r_def), parameter :: dz = 2000.0_r_def - integer(i_def), parameter :: nlayers = 3 + real(r_def), parameter :: tol = 1.0e-6_r_def + real(r_def), parameter :: dx = 6000.0_r_def + real(r_def), parameter :: dy = 1000.0_r_def + real(r_def), parameter :: dz = 2000.0_r_def + integer(i_def), parameter :: nlayers = 3 + + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi - integer :: i, j, k, cell - integer(i_def) :: err + integer :: cell integer(i_def) :: ndf_w0, undf_w0, ndf_w1, ndf_w2, ndf_pid, undf_pid integer(i_def) :: ncells, ncell_3d, nqp_h, nqp_v @@ -164,6 +139,10 @@ contains chi2, & chi3, & panel_id, & + geometry, & + topology, & + coord_system, & + scaled_radius, & ndf_w2, & basis_w2, & ndf_w1, & diff --git a/components/science/unit-test/kernel/fem/compute_derham_matrices_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/compute_derham_matrices_kernel_mod_test.pf index 1674c37e8..31e7f52af 100644 --- a/components/science/unit-test/kernel/fem/compute_derham_matrices_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_derham_matrices_kernel_mod_test.pf @@ -8,77 +8,50 @@ !> module compute_derham_matrices_kernel_mod_test - use constants_mod, only: i_def, r_def, imdi + use constants_mod, only: i_def, r_def, imdi, rmdi + use finite_element_config_mod, only: coord_system_xyz use funit implicit none private - public :: test_all - - @TestCase - type, extends(TestCase), public :: compute_derham_matrices_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_derham_matrices_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_derham_matrices_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_derham_matrices_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_all( this ) + subroutine test_all() - use sci_compute_derham_matrices_kernel_mod, only : compute_derham_matrices_code + use sci_compute_derham_matrices_kernel_mod, only: compute_derham_matrices_code implicit none - class(compute_derham_matrices_test_type), intent(inout) :: this - real(r_def), parameter :: tol = 1.0e-12_r_def ! Operators to test @@ -107,6 +80,11 @@ contains integer(i_def) :: i real(kind=r_def) :: answer + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + ! Set cell coordinates dx = 2.0_r_def dy = 3.0_r_def @@ -157,6 +135,9 @@ contains 1_i_def, broken_div, & x, y, z, & panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf, s, ds, & ndf, v, ds, & ndf, v, dv, & diff --git a/components/science/unit-test/kernel/fem/compute_div_operator_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/compute_div_operator_kernel_mod_test.pf index 95e85be0a..9f48a1452 100644 --- a/components/science/unit-test/kernel/fem/compute_div_operator_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_div_operator_kernel_mod_test.pf @@ -8,7 +8,9 @@ !> module compute_div_operator_kernel_mod_test - use constants_mod, only : i_def, r_def, imdi + use constants_mod, only : l_def, i_def, r_def, imdi, rmdi + use finite_element_config_mod, only: coord_system_xyz + use funit use get_unit_test_m3x3_dofmap_mod, & only : get_w0_m3x3_dofmap, & get_w3_m3x3_dofmap @@ -26,82 +28,59 @@ module compute_div_operator_kernel_mod_test get_gaussian_q3x3x3_quadrature_weights_z use get_unit_test_3x3x3_chi_mod, & only : get_w0_3x3x3_field - use funit implicit none private - public :: test_all - - @TestCase - type, extends(TestCase), public :: compute_div_operator_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_div_operator_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_div_operator_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=1_i_def, & - element_order_v=1_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_div_operator_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_all( this ) + subroutine test_all() use sci_compute_div_operator_kernel_mod, only : compute_div_operator_code implicit none - class(compute_div_operator_test_type), intent(inout) :: this - real(kind=r_def), parameter :: tol = 1.0e-6_r_def real(kind=r_def), parameter :: dx = 6000.0_r_def real(r_def), parameter :: dy = 1000.0_r_def real(r_def), parameter :: dz = 2000.0_r_def integer(i_def), parameter :: nlayers = 3 - integer :: i, j, k, cell + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + logical(l_def), parameter :: rehabilitate = .true. + + integer :: cell integer(i_def) :: ndf_w0, undf_w0, ndf_w2, ndf_w3, ndf_pid, undf_pid integer(i_def) :: nqp_h, nqp_v @@ -165,6 +144,11 @@ contains chi2, & chi3, & panel_id, & + geometry, & + topology, & + coord_system, & + scaled_radius, & + rehabilitate, & ndf_w3, & basis_w3, & ndf_w2, & diff --git a/components/science/unit-test/kernel/fem/compute_grad_operator_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/compute_grad_operator_kernel_mod_test.pf index 3ee9177c2..57c7b4ecc 100644 --- a/components/science/unit-test/kernel/fem/compute_grad_operator_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_grad_operator_kernel_mod_test.pf @@ -8,7 +8,9 @@ !> module compute_grad_operator_kernel_mod_test - use constants_mod, only : i_def, r_def, imdi + use constants_mod, only : i_def, r_def, imdi, rmdi + use finite_element_config_mod, only : coord_system_xyz + use funit use get_unit_test_m3x3_q3x3x3_sizes_mod, only : get_w0_m3x3_q3x3x3_size, & get_w1_m3x3_q3x3x3_size, & get_w3_m3x3_q3x3x3_size @@ -20,91 +22,64 @@ module compute_grad_operator_kernel_mod_test use get_unit_test_q3x3x3_quadrature_mod, only : get_gaussian_q3x3x3_quadrature_weights_xy, & get_gaussian_q3x3x3_quadrature_weights_z use get_unit_test_3x3x3_chi_mod, only : get_w0_3x3x3_field - use funit implicit none private - public :: test_all - - @TestCase - type, extends(TestCase), public :: compute_grad_operator_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_grad_operator_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_grad_operator_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_grad_operator_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - @Test - subroutine test_all( this ) + subroutine test_all() use sci_compute_grad_operator_kernel_mod, & only : compute_grad_operator_code implicit none - class(compute_grad_operator_test_type), intent(inout) :: this - real(kind=r_def), dimension(12,8) :: answer real(kind=r_def), parameter :: tol = 1.0e-6_r_def real(kind=r_def), parameter :: dx = 6000.0_r_def real(r_def), parameter :: dy = 1000.0_r_def real(r_def), parameter :: dz = 2000.0_r_def - integer(i_def) :: nlayers, ncells, ncell_3d, cell, icell + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + + integer(i_def) :: nlayers, ncells, ncell_3d, cell integer(i_def) :: ndf_w0, undf_w0, ndf_w1, undf_w1 integer(i_def) :: ndf_pid, undf_pid, ndf_wchi, undf_wchi integer(i_def) :: dim_space, dim_space_diff integer(i_def) :: nqp_h, nqp_v - integer(i_def) :: i, j, k - integer(i_def), allocatable :: map_wchi(:,:), map_pid(:,:) real(r_def), allocatable :: diff_basis_w0(:,:,:,:), & basis_w1(:,:,:,:), & @@ -167,6 +142,10 @@ contains chi2, & chi3, & panel_id, & + geometry, & + topology, & + coord_system, & + scaled_radius, & ndf_w1, & basis_w1, & ndf_w0, & diff --git a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w0_mod_test.pf b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w0_mod_test.pf index ddaaa794d..89795c1d5 100644 --- a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w0_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w0_mod_test.pf @@ -8,7 +8,7 @@ !> module compute_mass_matrix_kernel_w0_mod_test - use constants_mod, only : i_def, r_def, r_single, r_double, l_def + use constants_mod, only: i_def, r_def, imdi, rmdi, r_single, r_double, l_def use get_unit_test_m3x3_q3x3x3_sizes_mod, only : get_w0_m3x3_q3x3x3_size, & get_w3_m3x3_q3x3x3_size @@ -22,72 +22,25 @@ module compute_mass_matrix_kernel_w0_mod_test use get_unit_test_3x3x3_chi_mod, only : get_w0_3x3x3_field use funit + use finite_element_config_mod, only: coord_system_xyz + implicit none private public :: test_r_single, test_mixed_precision, test_r_double - @TestCase - type, extends(TestCase), public :: compute_mass_matrix_w0_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_r_single - procedure test_mixed_precision - procedure test_r_double - end type compute_mass_matrix_w0_test_type - contains - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) - - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - - implicit none - - class(compute_mass_matrix_w0_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape = cellshape_quadrilateral, & - coord_order = 0_i_def, & - coord_order_nonprime = 1_i_def, & - coord_space = coord_space_wchi, & - coord_system = coord_system_xyz, & - element_order_h =0_i_def, & - element_order_v =0_i_def, & - rehabilitate = .true. ) - - end subroutine setUp - - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) - - use config_loader_mod, only: final_configuration - - implicit none - - class(compute_mass_matrix_w0_test_type), intent(inout) :: this - - call final_configuration() - - end subroutine tearDown - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! TEST R_SINGLE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_r_single( this ) + subroutine test_r_single() use sci_compute_mass_matrix_kernel_w_scalar_mod, only : & compute_mass_matrix_w_scalar_code implicit none - class(compute_mass_matrix_w0_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def @@ -120,6 +73,10 @@ contains real(r_single), allocatable :: diff_basis_w0_r32(:,:,:,:) real(r_single), allocatable :: wh_r32(:), wv_r32(:) + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi ! Get pre-set test data (in i_def and r_def) call get_w0_m3x3_q3x3x3_size( ndf_w0, undf_w0, ncells, & @@ -187,6 +144,9 @@ contains chi3_r32, & panel_id_r32, & extended_mesh, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf_w0, & basis_w0_r32, & ndf_w0, & @@ -227,15 +187,13 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!! TEST MIXED PRECISION !!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_mixed_precision( this ) + subroutine test_mixed_precision() use sci_compute_mass_matrix_kernel_w_scalar_mod, only : & compute_mass_matrix_w_scalar_code implicit none - class(compute_mass_matrix_w0_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def @@ -268,6 +226,10 @@ contains real(r_double), allocatable :: diff_basis_w0_r64(:,:,:,:) real(r_double), allocatable :: wh_r64(:), wv_r64(:) + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi ! Get pre-set test data (in i_def and r_def) call get_w0_m3x3_q3x3x3_size( ndf_w0, undf_w0, ncells, & @@ -335,6 +297,9 @@ contains chi3_r64, & panel_id_r64, & extended_mesh, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf_w0, & basis_w0_r64, & ndf_w0, & @@ -375,15 +340,13 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! TEST R_DOUBLE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_r_double( this ) + subroutine test_r_double() use sci_compute_mass_matrix_kernel_w_scalar_mod, only : & compute_mass_matrix_w_scalar_code implicit none - class(compute_mass_matrix_w0_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def @@ -416,6 +379,10 @@ contains real(r_double), allocatable :: diff_basis_w0_r64(:,:,:,:) real(r_double), allocatable :: wh_r64(:), wv_r64(:) + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi ! Get pre-set test data (in i_def and r_def) call get_w0_m3x3_q3x3x3_size( ndf_w0, undf_w0, ncells, & @@ -483,6 +450,9 @@ contains chi3_r64, & panel_id_r64, & extended_mesh, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf_w0, & basis_w0_r64, & ndf_w0, & diff --git a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w1_mod_test.pf b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w1_mod_test.pf index 204233951..0e03791c7 100644 --- a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w1_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w1_mod_test.pf @@ -8,94 +8,73 @@ !> module compute_mass_matrix_kernel_w1_mod_test - use constants_mod, only : i_def, r_def, imdi - use get_unit_test_m3x3_dofmap_mod, & + use constants_mod, only : i_def, r_def, imdi, rmdi + use get_unit_test_m3x3_dofmap_mod, & only : get_w0_m3x3_dofmap, get_w3_m3x3_dofmap - use get_unit_test_m3x3_q3x3x3_sizes_mod, & + use get_unit_test_m3x3_q3x3x3_sizes_mod, & only : get_w0_m3x3_q3x3x3_size, get_w1_m3x3_q3x3x3_size, & get_w3_m3x3_q3x3x3_size - use get_unit_test_q3x3x3_basis_mod, & + use get_unit_test_q3x3x3_basis_mod, & only : get_w0_q3x3x3_diff_basis, get_w1_q3x3x3_basis, & get_w0_q3x3x3_basis - use get_unit_test_q3x3x3_quadrature_mod, & + use get_unit_test_q3x3x3_quadrature_mod, & only : get_gaussian_q3x3x3_quadrature_weights_xy, & get_gaussian_q3x3x3_quadrature_weights_z - use get_unit_test_3x3x3_chi_mod, & + use get_unit_test_3x3x3_chi_mod, & only : get_w0_3x3x3_field - use funit + use funit + + use finite_element_config_mod, only: coord_system_xyz implicit none private - public :: test_all - - @TestCase - type, extends(TestCase), public :: compute_mass_matrix_w1_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_mass_matrix_w1_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_mass_matrix_w1_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_mass_matrix_w1_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_all( this ) + subroutine test_all() use sci_compute_mass_matrix_kernel_w1_mod, only : compute_mass_matrix_w1_code implicit none - class(compute_mass_matrix_w1_test_type), intent(inout) :: this - real(r_def), parameter :: tol = 1.0e-6_r_def real(r_def), parameter :: dx = 6000.0_r_def real(r_def), parameter :: dy = 1000.0_r_def real(r_def), parameter :: dz = 2000.0_r_def + integer(i_def), parameter :: nlayers = 3 + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi integer :: i, j, k, cell @@ -159,6 +138,10 @@ contains chi2, & chi3, & panel_id, & + geometry, & + topology, & + coord_system, & + scaled_radius, & ndf_w1, & basis_w1, & ndf_w0, & diff --git a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2_mod_test.pf b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2_mod_test.pf index 9fe32c78b..f0ceaa119 100644 --- a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2_mod_test.pf @@ -8,7 +8,8 @@ !> module compute_mass_matrix_kernel_w2_mod_test - use constants_mod, only : i_def, r_def, imdi + use constants_mod, only: i_def, r_def, imdi, rmdi + use get_unit_test_m3x3_dofmap_mod, & only : get_w0_m3x3_dofmap, get_w3_m3x3_dofmap use get_unit_test_m3x3_q3x3x3_sizes_mod, & @@ -24,78 +25,57 @@ module compute_mass_matrix_kernel_w2_mod_test only : get_w0_3x3x3_field use funit + use finite_element_config_mod, only: coord_system_xyz + implicit none private - public :: test_all - - @TestCase - type, extends(TestCase), public :: compute_mass_matrix_w2_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_mass_matrix_w2_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_mass_matrix_w2_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_mass_matrix_w2_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_all( this ) + subroutine test_all() use sci_compute_mass_matrix_kernel_w2_mod, only : compute_mass_matrix_w2_code implicit none - class(compute_mass_matrix_w2_test_type), intent(inout) :: this - real(r_def), parameter :: tol = 1.0e-6_r_def real(r_def), parameter :: dx = 6000.0_r_def real(r_def), parameter :: dy = 1000.0_r_def real(r_def), parameter :: dz = 2000.0_r_def + integer(i_def), parameter :: nlayers = 3 + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi integer :: i, j, k, cell @@ -160,6 +140,10 @@ contains chi2, & chi3, & panel_id, & + geometry, & + topology, & + coord_system, & + scaled_radius, & ndf_w2, & basis_w2, & ndf_w0, & diff --git a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2b_mod_test.pf b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2b_mod_test.pf index fe9c58026..043fdfdf2 100644 --- a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2b_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w2b_mod_test.pf @@ -7,7 +7,8 @@ !> module compute_mass_matrix_kernel_w2b_mod_test - use constants_mod, only : i_def, r_def + use constants_mod, only: i_def, r_def, imdi, rmdi + use get_unit_test_m3x3_q3x3x3_sizes_mod, only : get_w0_m3x3_q3x3x3_size, & get_w2broken_m3x3_q3x3x3_size, & get_w3_m3x3_q3x3x3_size @@ -21,75 +22,33 @@ module compute_mass_matrix_kernel_w2b_mod_test use get_unit_test_3x3x3_chi_mod, only : get_w0_3x3x3_field use funit + use finite_element_config_mod, only: coord_system_xyz + implicit none private public :: test_all - @TestCase - type, extends(TestCase), public :: compute_mass_matrix_w2b_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type compute_mass_matrix_w2b_test_type - contains - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) - - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - - implicit none - - class(compute_mass_matrix_w2b_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - - end subroutine setUp - - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) - - use config_loader_mod, only: final_configuration - - implicit none - - class(compute_mass_matrix_w2b_test_type), intent(inout) :: this - - call final_configuration() - - end subroutine tearDown - - !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - @Test - subroutine test_all( this ) + subroutine test_all() use sci_compute_mass_matrix_kernel_w2_mod, only : compute_mass_matrix_w2_code implicit none - class(compute_mass_matrix_w2b_test_type), intent(inout) :: this - real(kind=r_def) :: answer(6,6) real(r_def), parameter :: tol = 1.0e-6_r_def real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + integer(i_def) :: nlayers, ncells, ncell_3d, cell, icell integer(i_def) :: ndf_w0, undf_w0, ndf_w2broken, undf_w2broken, ndf_w3, undf_w3 integer(i_def) :: dim_space, dim_space_diff @@ -151,6 +110,10 @@ contains chi2, & chi3, & panel_id, & + geometry, & + topology, & + coord_system, & + scaled_radius, & ndf_w2broken, & basis_w2broken, & ndf_w0, & diff --git a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w3_mod_test.pf b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w3_mod_test.pf index 8a9cedd13..9bf127e7b 100644 --- a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w3_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_w3_mod_test.pf @@ -8,7 +8,7 @@ !> module compute_mass_matrix_kernel_w3_mod_test - use constants_mod, only : i_def, r_def, r_single, r_double, imdi + use constants_mod, only: i_def, r_def, l_def, r_single, r_double, imdi, rmdi use get_unit_test_m3x3_q3x3x3_sizes_mod, only : get_w0_m3x3_q3x3x3_size, & get_w3_m3x3_q3x3x3_size @@ -22,81 +22,61 @@ module compute_mass_matrix_kernel_w3_mod_test use get_unit_test_3x3x3_chi_mod, only : get_w0_3x3x3_field use funit + use finite_element_config_mod, only: coord_system_xyz + implicit none private + public :: set_up, tear_down public :: test_r_single, test_r_double, test_mixed_precision - @TestCase - type, extends(TestCase), public :: compute_mass_matrix_w3_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_r_single - procedure test_r_double - procedure test_mixed_precision - end type compute_mass_matrix_w3_test_type - contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_mass_matrix_w3_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape = cellshape_quadrilateral, & - coord_order = 0_i_def, & - coord_order_nonprime = 1_i_def, & - coord_space = coord_space_wchi, & - coord_system = coord_system_xyz, & - element_order_h = 0_i_def, & - element_order_v = 0_i_def, & - rehabilitate = .true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_mass_matrix_w3_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!! TEST R_SINGLE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_r_single( this ) + subroutine test_r_single() use sci_compute_mass_matrix_kernel_w3_mod, only : compute_mass_matrix_w3_code implicit none - class(compute_mass_matrix_w3_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def integer(i_def), parameter :: nlayers = 3_i_def real(r_single), parameter :: tol = 1.0e-6_r_single + + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + logical(l_def), parameter :: rehabilitate = .true. + real(r_single) :: answer integer(i_def) :: ncells, ncell_3d, cell @@ -198,6 +178,10 @@ contains chi2_r32, & chi3_r32, & panel_id_r32, & + geometry, topology, & + coord_system, & + scaled_radius, & + rehabilitate, & ndf_w3, & basis_w3_r32, & ndf_w0, & @@ -234,19 +218,24 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!! TEST R_DOUBLE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_r_double( this ) + subroutine test_r_double() use sci_compute_mass_matrix_kernel_w3_mod, only : compute_mass_matrix_w3_code implicit none - class(compute_mass_matrix_w3_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def integer(i_def), parameter :: nlayers = 3_i_def real(r_double), parameter :: tol = 1.0e-6_r_double + + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + logical(l_def), parameter :: rehabilitate = .true. + real(r_double) :: answer integer(i_def) :: ncells, ncell_3d, cell @@ -348,6 +337,10 @@ contains chi2_r64, & chi3_r64, & panel_id_r64, & + geometry, topology, & + coord_system, & + scaled_radius, & + rehabilitate, & ndf_w3, & basis_w3_r64, & ndf_w0, & @@ -385,19 +378,24 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!! TEST MIXED PRECISION !!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_mixed_precision( this ) + subroutine test_mixed_precision() use sci_compute_mass_matrix_kernel_w3_mod, only : compute_mass_matrix_w3_code implicit none - class(compute_mass_matrix_w3_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def integer(i_def), parameter :: nlayers = 3_i_def real(r_single), parameter :: tol = 1.0e-6_r_single + + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + logical(l_def), parameter :: rehabilitate = .true. + real(r_single) :: answer integer(i_def) :: ncells, ncell_3d, cell @@ -499,6 +497,10 @@ contains chi2_r64, & chi3_r64, & panel_id_r64, & + geometry, topology, & + coord_system, & + scaled_radius, & + rehabilitate, & ndf_w3, & basis_w3_r64, & ndf_w0, & diff --git a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_wtheta_mod_test.pf b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_wtheta_mod_test.pf index 43334a4fe..b09876167 100644 --- a/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_wtheta_mod_test.pf +++ b/components/science/unit-test/kernel/fem/compute_mass_matrix_kernel_wtheta_mod_test.pf @@ -8,7 +8,8 @@ !> module compute_mass_matrix_kernel_wtheta_mod_test - use constants_mod, only : i_def, r_def, r_single, r_double, l_def, imdi + use constants_mod, only: i_def, r_def, r_single, r_double, & + l_def, imdi, rmdi use get_unit_test_m3x3_q3x3x3_sizes_mod, & only : get_w0_m3x3_q3x3x3_size, & @@ -28,77 +29,50 @@ module compute_mass_matrix_kernel_wtheta_mod_test only : get_w0_3x3x3_field use funit + use finite_element_config_mod, only: coord_system_xyz + implicit none private + public set_up, tear_down public :: test_r_single, test_mixed_precision, test_r_double - @TestCase - type, extends(TestCase), public :: compute_mass_matrix_wtheta_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_r_single - procedure test_mixed_precision - procedure test_r_double - end type compute_mass_matrix_wtheta_test_type - contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(compute_mass_matrix_wtheta_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape = cellshape_quadrilateral, & - coord_order = 0_i_def, & - coord_order_nonprime = 1_i_def, & - coord_space = coord_space_wchi, & - coord_system = coord_system_xyz, & - element_order_h = 0_i_def, & - element_order_v = 0_i_def, & - rehabilitate = .true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + @after + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(compute_mass_matrix_wtheta_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! TEST R_SINGLE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_r_single( this ) + subroutine test_r_single() use sci_compute_mass_matrix_kernel_w_scalar_mod, only : & compute_mass_matrix_w_scalar_code implicit none - class(compute_mass_matrix_wtheta_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def @@ -118,6 +92,10 @@ contains real(r_single), allocatable :: local_stencil_r32(:,:,:) logical(l_def), parameter :: extended_mesh = .false. + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi ! Canned test data is kept as r_def real(r_def), allocatable :: chi1(:), chi2(:), chi3(:) @@ -211,6 +189,9 @@ contains chi3_r32, & panel_id_r32, & extended_mesh, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf_wtheta, & basis_wtheta_r32, & ndf_w0, & @@ -255,15 +236,13 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!! TEST MIXED PRECISION !!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_mixed_precision( this ) + subroutine test_mixed_precision() use sci_compute_mass_matrix_kernel_w_scalar_mod, only : & compute_mass_matrix_w_scalar_code implicit none - class(compute_mass_matrix_wtheta_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def @@ -283,6 +262,10 @@ contains real(r_single), allocatable :: local_stencil_r32(:,:,:) logical(l_def), parameter :: extended_mesh = .false. + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi ! Canned test data is kept as r_def real(r_def), allocatable :: chi1(:), chi2(:), chi3(:) @@ -376,6 +359,9 @@ contains chi3_r64, & panel_id_r64, & extended_mesh, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf_wtheta, & basis_wtheta_r64, & ndf_w0, & @@ -420,15 +406,13 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! TEST R_DOUBLE !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! @Test - subroutine test_r_double( this ) + subroutine test_r_double() use sci_compute_mass_matrix_kernel_w_scalar_mod, only : & compute_mass_matrix_w_scalar_code implicit none - class(compute_mass_matrix_wtheta_test_type), intent(inout) :: this - real(r_def), parameter :: dx = 6000.0_r_def, & dy = 1000.0_r_def, & dz = 2000.0_r_def @@ -448,6 +432,10 @@ contains real(r_double), allocatable :: local_stencil_r64(:,:,:) logical(l_def), parameter :: extended_mesh = .false. + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi ! Canned test data is kept as r_def real(r_def), allocatable :: chi1(:), chi2(:), chi3(:) @@ -541,6 +529,9 @@ contains chi3_r64, & panel_id_r64, & extended_mesh, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf_wtheta, & basis_wtheta_r64, & ndf_w0, & diff --git a/components/science/unit-test/kernel/fem/gp_rhs_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/gp_rhs_kernel_mod_test.pf index e3e582d9f..8c880313b 100644 --- a/components/science/unit-test/kernel/fem/gp_rhs_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/gp_rhs_kernel_mod_test.pf @@ -7,11 +7,12 @@ !------------------------------------------------------------------------------- module gp_rhs_kernel_mod_test - use constants_mod, only : i_def, r_def, imdi + use constants_mod, only: i_def, r_def, imdi, rmdi + use finite_element_config_mod, only: coord_system_xyz use funit - use quadrature_xyoz_mod, only: quadrature_xyoz_type, & - quadrature_xyoz_proxy_type - use quadrature_rule_gaussian_mod, only: quadrature_rule_gaussian_type + use quadrature_xyoz_mod, only: quadrature_xyoz_type, & + quadrature_xyoz_proxy_type + use quadrature_rule_gaussian_mod, only: quadrature_rule_gaussian_type implicit none @@ -40,33 +41,21 @@ module gp_rhs_kernel_mod_test contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + @before subroutine setUp( this ) - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none class(gp_rhs_test_type), intent(inout) :: this type( quadrature_xyoz_type ) :: qr - type(quadrature_xyoz_proxy_type) :: qr_proxy + type( quadrature_xyoz_proxy_type ) :: qr_proxy type( quadrature_rule_gaussian_type ) :: quadrature_rule integer(i_def) :: nqp_h, nqp_v - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) qr = quadrature_xyoz_type(3, quadrature_rule) @@ -75,16 +64,16 @@ contains nqp_h=qr_proxy%np_xy nqp_v=qr_proxy%np_z - allocate( this%basis_f (1,ndf_f,nqp_h,nqp_v), & - this%basis (1,ndf_w0,nqp_h,nqp_v), & + allocate( this%basis_f(1,ndf_f,nqp_h,nqp_v), & + this%basis(1,ndf_w0,nqp_h,nqp_v), & this%diff_basis(3,ndf_w0,nqp_h,nqp_v) ) end subroutine setUp !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + @after subroutine tearDown( this ) - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none @@ -93,7 +82,6 @@ contains deallocate( this%basis, this%basis_f, this%diff_basis ) - call final_configuration() call final_chi_transforms() end subroutine tearDown @@ -124,6 +112,11 @@ contains real(kind=r_def) :: fx, fy, fz integer :: df, dim, qp1, qp2 + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + qr = quadrature_xyoz_type(3, quadrature_rule) qr_proxy = qr%get_quadrature_proxy() @@ -150,6 +143,7 @@ contains fx = (1.0_r_def - x_vert(df,1)) + (-1.0_r_def)**(int(x_vert(df,1))+1)*xqp_h(qp1,1) fy = (1.0_r_def - x_vert(df,2)) + (-1.0_r_def)**(int(x_vert(df,2))+1)*xqp_h(qp1,2) fz = (1.0_r_def - x_vert(df,3)) + (-1.0_r_def)**(int(x_vert(df,3))+1)*xqp_v(qp2) + this%diff_basis(1,df,qp1,qp2) = fy*fz this%diff_basis(2,df,qp1,qp2) = fx*fz this%diff_basis(3,df,qp1,qp2) = fx*fy @@ -172,15 +166,16 @@ contains num_layers=1 rhs(:) = 0.0_r_def - call gp_rhs_code(num_layers, & - rhs, f_data, & - chi(:,1), chi(:,2), chi(:,3), & - panel_id, & - ndf_w0, undf, map_w0, this%basis, & - ndf_f, undf_f, map_f, this%basis_f, & - ndf_w0, undf, map_w0, & - this%basis, this%diff_basis, & - ndf_pid, undf_pid, map_pid, & + call gp_rhs_code(num_layers, & + rhs, f_data, & + chi(:,1), chi(:,2), chi(:,3), & + panel_id, geometry, topology, & + coord_system, scaled_radius, & + ndf_w0, undf, map_w0, this%basis, & + ndf_f, undf_f, map_f, this%basis_f, & + ndf_w0, undf, map_w0, & + this%basis, this%diff_basis, & + ndf_pid, undf_pid, map_pid, & nqp_h, nqp_v, wh, wv ) ! Answer for hardwired jacobian diff --git a/components/science/unit-test/kernel/fem/gp_vector_rhs_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/gp_vector_rhs_kernel_mod_test.pf index 88cd570cb..55eb77f53 100644 --- a/components/science/unit-test/kernel/fem/gp_vector_rhs_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/gp_vector_rhs_kernel_mod_test.pf @@ -7,11 +7,15 @@ !------------------------------------------------------------------------------- module gp_vector_rhs_kernel_mod_test - use constants_mod, only : i_def, r_def + use constants_mod, only: i_def, r_def, rmdi + use finite_element_config_mod, only: coord_system_xyz use funit - use quadrature_xyoz_mod, only: quadrature_xyoz_type, & - quadrature_xyoz_proxy_type - use quadrature_rule_gaussian_mod, only: quadrature_rule_gaussian_type + use quadrature_xyoz_mod, only: quadrature_xyoz_type, & + quadrature_xyoz_proxy_type + use quadrature_rule_gaussian_mod, only: quadrature_rule_gaussian_type + + use sci_mesh_enums_mod, only: geometry_planar, & + topology_fully_periodic implicit none @@ -36,43 +40,21 @@ module gp_vector_rhs_kernel_mod_test contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + @before subroutine setUp( this ) - use base_mesh_config_mod, only : geometry_planar, & - topology_fully_periodic - use sci_chi_transform_mod, only : init_chi_transforms - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config, & - feign_base_mesh_config + use sci_chi_transform_mod, only: init_chi_transforms implicit none class(gp_vector_rhs_test_type), intent(inout) :: this type( quadrature_xyoz_type ) :: qr - type(quadrature_xyoz_proxy_type) :: qr_proxy + type( quadrature_xyoz_proxy_type ) :: qr_proxy type( quadrature_rule_gaussian_type ) :: quadrature_rule integer(i_def) :: nqp_h, nqp_v - call feign_base_mesh_config( file_prefix='foo', & - prime_mesh_name='unit_test', & - geometry=geometry_planar, & - prepartitioned=.false., & - topology=topology_fully_periodic, & - fplane=.false., f_lat_deg=0.0_r_def ) - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(geometry_planar, topology_fully_periodic) qr = quadrature_xyoz_type(3, quadrature_rule) @@ -88,9 +70,9 @@ contains end subroutine setUp !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! + @after subroutine tearDown( this ) - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none @@ -99,7 +81,6 @@ contains deallocate( this%basis_f, this%basis, this%diff_basis ) - call final_configuration() call final_chi_transforms() end subroutine tearDown @@ -129,6 +110,10 @@ contains real(kind=r_def) :: fx, fy, fz, dx, dy, dz, t_vec(3) integer :: df, dim, qp1, qp2 + integer(i_def), parameter :: geometry = geometry_planar + integer(i_def), parameter :: topology = topology_fully_periodic + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi qr = quadrature_xyoz_type(3, quadrature_rule) @@ -193,6 +178,10 @@ contains chi(:,1), chi(:,2), chi(:,3), & panel_id, & f_data, & + geometry, & + topology, & + coord_system, & + scaled_radius, & ndf_w0, undf, map_w0, this%basis, & ndf_f, undf_f,map_f, this%basis_f, & ndf_w0, undf, map_w0, this%basis, & diff --git a/components/science/unit-test/kernel/fem/mg_derham_mat_kernel_mod_test.pf b/components/science/unit-test/kernel/fem/mg_derham_mat_kernel_mod_test.pf index 3ac46475e..49dc19a8b 100644 --- a/components/science/unit-test/kernel/fem/mg_derham_mat_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/fem/mg_derham_mat_kernel_mod_test.pf @@ -8,75 +8,48 @@ !> module mg_derham_mat_kernel_mod_test - use constants_mod, only : i_def, r_def, imdi + use constants_mod, only : i_def, r_def, imdi, rmdi + use finite_element_config_mod, only: coord_system_xyz + use funit implicit none private - public :: test_all - - @TestCase - type, extends(TestCase), public :: mg_derham_mat_test_type - private - contains - procedure setUp - procedure tearDown - procedure test_all - end type mg_derham_mat_test_type + public :: set_up, tear_down, test_all contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine setUp( this ) + @before + subroutine set_up() - use finite_element_config_mod, only : cellshape_quadrilateral, & - coord_system_xyz, coord_space_wchi - use feign_config_mod, only : feign_finite_element_config - use sci_chi_transform_mod, only : init_chi_transforms + use sci_chi_transform_mod, only: init_chi_transforms implicit none - class(mg_derham_mat_test_type), intent(inout) :: this - - call feign_finite_element_config( & - cellshape=cellshape_quadrilateral, & - coord_order=0_i_def, & - coord_order_nonprime=1_i_def, & - coord_space=coord_space_wchi, & - coord_system=coord_system_xyz, & - element_order_h=0_i_def, & - element_order_v=0_i_def, & - rehabilitate=.true. ) - call init_chi_transforms(imdi, imdi) - end subroutine setUp + end subroutine set_up !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! - subroutine tearDown( this ) + subroutine tear_down() - use config_loader_mod, only: final_configuration use sci_chi_transform_mod, only: final_chi_transforms implicit none - class(mg_derham_mat_test_type), intent(inout) :: this - - call final_configuration() call final_chi_transforms() - end subroutine tearDown + end subroutine tear_down @Test - subroutine test_all( this ) + subroutine test_all() use sci_mg_derham_mat_kernel_mod, only : mg_derham_mat_code implicit none - class(mg_derham_mat_test_type), intent(inout) :: this - real(r_def), parameter :: tol = 1.0e-12_r_def ! Operators to test @@ -105,6 +78,11 @@ contains integer(i_def) :: i real(kind=r_def) :: answer + integer(i_def), parameter :: geometry = imdi + integer(i_def), parameter :: topology = imdi + integer(i_def), parameter :: coord_system = coord_system_xyz + real(r_def), parameter :: scaled_radius = rmdi + ! Set cell coordinates dx = 2.0_r_def dy = 3.0_r_def @@ -149,6 +127,9 @@ contains 1_i_def, div, & x, y, z, & panel_id, & + geometry, topology, & + coord_system, & + scaled_radius, & ndf, v, dv, & ndf, s, & ndf, s, & diff --git a/components/science/unit-test/kernel/geometry/chi_transform_mod_test.pf b/components/science/unit-test/kernel/geometry/chi_transform_mod_test.pf index 26f584216..23d0577a9 100644 --- a/components/science/unit-test/kernel/geometry/chi_transform_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/chi_transform_mod_test.pf @@ -10,14 +10,13 @@ module chi_transform_mod_test use constants_mod, only : i_def, r_def, str_long, PI, rmdi - use base_mesh_config_mod, only: geometry_spherical, & - geometry_planar, & - topology_fully_periodic, & - topology_non_periodic - use finite_element_config_mod, only: coord_system_native, & coord_system_xyz + use sci_mesh_enums_mod, only: geometry_spherical, & + geometry_planar, & + topology_fully_periodic, & + topology_non_periodic use funit implicit none diff --git a/components/science/unit-test/kernel/geometry/coordinate_jacobian_alphabetaz_mod_test.pf b/components/science/unit-test/kernel/geometry/coordinate_jacobian_alphabetaz_mod_test.pf index ae1aa557a..71d26f2b4 100644 --- a/components/science/unit-test/kernel/geometry/coordinate_jacobian_alphabetaz_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/coordinate_jacobian_alphabetaz_mod_test.pf @@ -10,9 +10,9 @@ module coordinate_jacobian_alphabetaz_mod_test use funit use constants_mod, only: r_def, i_def - use base_mesh_config_mod, only: geometry_spherical, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_native + use sci_mesh_enums_mod, only: geometry_spherical, & + topology_fully_periodic implicit none diff --git a/components/science/unit-test/kernel/geometry/coordinate_jacobian_lonlatz_mod_test.pf b/components/science/unit-test/kernel/geometry/coordinate_jacobian_lonlatz_mod_test.pf index e4a2c021f..c9ba5eb01 100644 --- a/components/science/unit-test/kernel/geometry/coordinate_jacobian_lonlatz_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/coordinate_jacobian_lonlatz_mod_test.pf @@ -9,11 +9,9 @@ module coordinate_jacobian_lonlatz_mod_test use funit use constants_mod, only : r_def, i_def, PI - - use base_mesh_config_mod, only: geometry_spherical, & - topology_non_periodic use finite_element_config_mod, only: coord_system_native - + use sci_mesh_enums_mod, only: geometry_spherical, & + topology_non_periodic implicit none public :: set_up, tear_down, test_all diff --git a/components/science/unit-test/kernel/geometry/coordinate_jacobian_stretched_mod_test.pf b/components/science/unit-test/kernel/geometry/coordinate_jacobian_stretched_mod_test.pf index abcec08a4..2b697b6dc 100644 --- a/components/science/unit-test/kernel/geometry/coordinate_jacobian_stretched_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/coordinate_jacobian_stretched_mod_test.pf @@ -21,8 +21,6 @@ contains @before subroutine set_up() - use base_mesh_config_mod, only : geometry_spherical, & - topology_non_periodic use extrusion_config_mod, only : method_uniform, & stretching_method_linear use finite_element_config_mod, only : cellshape_quadrilateral, & @@ -32,6 +30,8 @@ contains feign_extrusion_config, & feign_finite_element_config, & feign_planet_config + use sci_mesh_enums_mod, only : geometry_spherical, & + topology_non_periodic implicit none diff --git a/components/science/unit-test/kernel/geometry/coordinate_jacobian_xyz_mod_test.pf b/components/science/unit-test/kernel/geometry/coordinate_jacobian_xyz_mod_test.pf index bfc90c002..7c9c58597 100644 --- a/components/science/unit-test/kernel/geometry/coordinate_jacobian_xyz_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/coordinate_jacobian_xyz_mod_test.pf @@ -64,7 +64,7 @@ contains pointwise_coordinate_jacobian_inverse use finite_element_config_mod, only: coord_system_xyz - use base_mesh_config_mod, only: geometry_planar, & + use sci_mesh_enums_mod, only: geometry_planar, & topology_non_periodic implicit none diff --git a/components/science/unit-test/kernel/geometry/get_dz_w3_kernel_mod_test.pf b/components/science/unit-test/kernel/geometry/get_dz_w3_kernel_mod_test.pf index f3644f666..d4961e8a2 100644 --- a/components/science/unit-test/kernel/geometry/get_dz_w3_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/get_dz_w3_kernel_mod_test.pf @@ -35,12 +35,12 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine setUp( this ) - use base_mesh_config_mod, only : geometry_planar, & - topology_fully_periodic use finite_element_config_mod, only : cellshape_quadrilateral, & coord_system_xyz, coord_space_wchi use feign_config_mod, only : feign_finite_element_config, & feign_base_mesh_config + use sci_mesh_enums_mod, only : geometry_planar, & + topology_fully_periodic implicit none diff --git a/components/science/unit-test/kernel/geometry/height_continuous_kernel_mod_test.pf b/components/science/unit-test/kernel/geometry/height_continuous_kernel_mod_test.pf index b3e3d31e6..c19adf016 100644 --- a/components/science/unit-test/kernel/geometry/height_continuous_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/height_continuous_kernel_mod_test.pf @@ -7,8 +7,7 @@ !! This test is parametrised to test both planar and spherical geometries. module height_continuous_kernel_mod_test - use base_mesh_config_mod, only : geometry_planar, & - geometry_spherical + use constants_mod, only : i_def, r_def, str_long use funit use get_unit_test_w2nodal_basis_mod, only : get_w0_w2nodal_basis @@ -16,6 +15,8 @@ module height_continuous_kernel_mod_test get_w2_m3x3_q3x3x3_size use get_unit_test_m3x3_dofmap_mod, only : get_w0_m3x3_dofmap, & get_w2_m3x3_dofmap + use sci_mesh_enums_mod, only : geometry_planar, & + geometry_spherical implicit none diff --git a/components/science/unit-test/kernel/geometry/height_discontinuous_kernel_mod_test.pf b/components/science/unit-test/kernel/geometry/height_discontinuous_kernel_mod_test.pf index 400f083eb..8962268bd 100644 --- a/components/science/unit-test/kernel/geometry/height_discontinuous_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/height_discontinuous_kernel_mod_test.pf @@ -38,9 +38,9 @@ contains @Test subroutine test_all( this ) - use base_mesh_config_mod, only : geometry_planar use finite_element_config_mod, only : coord_system_xyz use sci_height_discontinuous_kernel_mod, only : height_discontinuous_code + use sci_mesh_enums_mod, only : geometry_planar implicit none diff --git a/components/science/unit-test/kernel/geometry/native_jacobian_alphabetaz_mod_test.pf b/components/science/unit-test/kernel/geometry/native_jacobian_alphabetaz_mod_test.pf index 162bc1642..f1d67375a 100644 --- a/components/science/unit-test/kernel/geometry/native_jacobian_alphabetaz_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/native_jacobian_alphabetaz_mod_test.pf @@ -9,9 +9,9 @@ module native_jacobian_alphabetaz_mod_test use funit use constants_mod, only: r_def, i_def - use base_mesh_config_mod, only: geometry_spherical, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_native + use sci_mesh_enums_mod, only: geometry_spherical, & + topology_fully_periodic implicit none diff --git a/components/science/unit-test/kernel/geometry/native_jacobian_lonlatz_mod_test.pf b/components/science/unit-test/kernel/geometry/native_jacobian_lonlatz_mod_test.pf index 53ed27607..f8c1707b0 100644 --- a/components/science/unit-test/kernel/geometry/native_jacobian_lonlatz_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/native_jacobian_lonlatz_mod_test.pf @@ -9,11 +9,9 @@ module native_jacobian_lonlatz_mod_test use funit use constants_mod, only : r_def, i_def, PI - - use base_mesh_config_mod, only: geometry_spherical, & - topology_non_periodic use finite_element_config_mod, only: coord_system_native - + use sci_mesh_enums_mod, only: geometry_spherical, & + topology_non_periodic implicit none public :: set_up, tear_down, test_all diff --git a/components/science/unit-test/kernel/geometry/native_jacobian_stretched_mod_test.pf b/components/science/unit-test/kernel/geometry/native_jacobian_stretched_mod_test.pf index f90c9b81b..0afaab188 100644 --- a/components/science/unit-test/kernel/geometry/native_jacobian_stretched_mod_test.pf +++ b/components/science/unit-test/kernel/geometry/native_jacobian_stretched_mod_test.pf @@ -21,8 +21,6 @@ contains @before subroutine set_up() - use base_mesh_config_mod, only : geometry_spherical, & - topology_non_periodic use extrusion_config_mod, only : method_uniform, & stretching_method_linear use finite_element_config_mod, only : cellshape_quadrilateral, & @@ -31,6 +29,8 @@ contains feign_extrusion_config, & feign_finite_element_config, & feign_planet_config + use sci_mesh_enums_mod, only : geometry_spherical, & + topology_non_periodic implicit none diff --git a/components/science/unit-test/kernel/inter_function_space/compute_map_u_operators_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_function_space/compute_map_u_operators_kernel_mod_test.pf index 5686d794b..45915f341 100644 --- a/components/science/unit-test/kernel/inter_function_space/compute_map_u_operators_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_function_space/compute_map_u_operators_kernel_mod_test.pf @@ -7,14 +7,10 @@ module compute_map_u_operators_kernel_mod_test use constants_mod, only : i_def, pi, r_def - - use base_mesh_config_mod, only: geometry_spherical, & - topology_non_periodic use finite_element_config_mod, only: coord_system_native - use funit - - + use sci_mesh_enums_mod, only: geometry_spherical, & + topology_non_periodic implicit none diff --git a/components/science/unit-test/kernel/inter_function_space/compute_sample_u_ops_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_function_space/compute_sample_u_ops_kernel_mod_test.pf index 53fe6db22..cbf9025ad 100644 --- a/components/science/unit-test/kernel/inter_function_space/compute_sample_u_ops_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_function_space/compute_sample_u_ops_kernel_mod_test.pf @@ -6,16 +6,12 @@ module compute_sample_u_ops_kernel_mod_test - use base_mesh_config_mod, only : geometry_spherical, & - topology_non_periodic use constants_mod, only : i_def, r_def use reference_element_mod, only : W, S, E, N, B, T - - use base_mesh_config_mod, only: geometry_spherical, & - topology_non_periodic use finite_element_config_mod, only: coord_system_native - use funit + use sci_mesh_enums_mod, only : geometry_spherical, & + topology_non_periodic implicit none diff --git a/components/science/unit-test/kernel/inter_function_space/convert_phys_to_hdiv_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_function_space/convert_phys_to_hdiv_kernel_mod_test.pf index 571181d93..5acf90c70 100644 --- a/components/science/unit-test/kernel/inter_function_space/convert_phys_to_hdiv_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_function_space/convert_phys_to_hdiv_kernel_mod_test.pf @@ -8,12 +8,10 @@ module convert_phys_to_hdiv_kernel_mod_test use constants_mod, only: i_def, r_def, rmdi - - use base_mesh_config_mod, only: geometry_planar, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_xyz - use funit + use sci_mesh_enums_mod, only: geometry_planar, & + topology_fully_periodic implicit none diff --git a/components/science/unit-test/kernel/inter_function_space/project_ws_to_w1_operator_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_function_space/project_ws_to_w1_operator_kernel_mod_test.pf index 5f1d0f284..2ee661987 100644 --- a/components/science/unit-test/kernel/inter_function_space/project_ws_to_w1_operator_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_function_space/project_ws_to_w1_operator_kernel_mod_test.pf @@ -23,9 +23,9 @@ module project_ws_to_w1_operator_kernel_mod_test get_w3_m3x3_dofmap use get_unit_test_3x3x3_chi_mod, only : get_w0_3x3x3_field - use base_mesh_config_mod, only: geometry_planar, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_xyz + use sci_mesh_enums_mod, only: geometry_planar, & + topology_fully_periodic implicit none diff --git a/components/science/unit-test/kernel/inter_function_space/sample_w3_to_wtheta_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_function_space/sample_w3_to_wtheta_kernel_mod_test.pf index 9f87b1d06..d10a0577b 100644 --- a/components/science/unit-test/kernel/inter_function_space/sample_w3_to_wtheta_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_function_space/sample_w3_to_wtheta_kernel_mod_test.pf @@ -34,10 +34,9 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine setUp( this ) - use base_mesh_config_mod, only : geometry_planar, & - topology_fully_periodic use feign_config_mod, only : feign_base_mesh_config - + use sci_mesh_enums_mod, only : geometry_planar, & + topology_fully_periodic implicit none class(sample_w3_to_wtheta_test_type), intent(inout) :: this diff --git a/components/science/unit-test/kernel/inter_function_space/sample_wtheta_to_w3_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_function_space/sample_wtheta_to_w3_kernel_mod_test.pf index 1f473d77a..d6c506baa 100644 --- a/components/science/unit-test/kernel/inter_function_space/sample_wtheta_to_w3_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_function_space/sample_wtheta_to_w3_kernel_mod_test.pf @@ -32,10 +32,9 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine setUp( this ) - use base_mesh_config_mod, only : geometry_planar, & - topology_fully_periodic use feign_config_mod, only : feign_base_mesh_config - + use sci_mesh_enums_mod, only : geometry_planar, & + topology_fully_periodic implicit none class(sample_wtheta_to_w3_test_type), intent(inout) :: this diff --git a/components/science/unit-test/kernel/inter_function_space/w3_to_w2_displacement_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_function_space/w3_to_w2_displacement_kernel_mod_test.pf index 6d01ea321..6682554f6 100644 --- a/components/science/unit-test/kernel/inter_function_space/w3_to_w2_displacement_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_function_space/w3_to_w2_displacement_kernel_mod_test.pf @@ -9,13 +9,10 @@ module w3_to_w2_displacement_kernel_mod_test use constants_mod, only: i_def, r_def, PI, l_def use reference_element_mod, only: S, E, N, W - - use base_mesh_config_mod, only: geometry_spherical, & - topology_fully_periodic use finite_element_config_mod, only: coord_system_native - use funit - + use sci_mesh_enums_mod, only: geometry_spherical, & + topology_fully_periodic implicit none diff --git a/components/science/unit-test/kernel/inter_mesh/inject_sh_w3_to_wt_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_mesh/inject_sh_w3_to_wt_kernel_mod_test.pf index 24cad0316..41563d793 100644 --- a/components/science/unit-test/kernel/inter_mesh/inject_sh_w3_to_wt_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_mesh/inject_sh_w3_to_wt_kernel_mod_test.pf @@ -32,7 +32,7 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine setUp( this ) - use base_mesh_config_mod, only : geometry_planar, & + use sci_mesh_enums_mod, only : geometry_planar, & topology_fully_periodic use feign_config_mod, only : feign_base_mesh_config diff --git a/components/science/unit-test/kernel/inter_mesh/inject_wt_to_sh_w3_kernel_mod_test.pf b/components/science/unit-test/kernel/inter_mesh/inject_wt_to_sh_w3_kernel_mod_test.pf index f5467afc2..903435156 100644 --- a/components/science/unit-test/kernel/inter_mesh/inject_wt_to_sh_w3_kernel_mod_test.pf +++ b/components/science/unit-test/kernel/inter_mesh/inject_wt_to_sh_w3_kernel_mod_test.pf @@ -32,9 +32,9 @@ contains !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!! subroutine setUp( this ) - use base_mesh_config_mod, only : geometry_planar, & - topology_fully_periodic use feign_config_mod, only : feign_base_mesh_config + use sci_mesh_enums_mod, only : geometry_planar, & + topology_fully_periodic implicit none