diff --git a/doc/pages/user-guide/simcomps.md b/doc/pages/user-guide/simcomps.md index 58372cdd95f..056214128d5 100644 --- a/doc/pages/user-guide/simcomps.md +++ b/doc/pages/user-guide/simcomps.md @@ -32,6 +32,7 @@ in Neko. The list will be updated as new simcomps are added. - Computation of the weak gradient of a field \ref simcomp_weak_grad - User defined components \ref user-file_simcomps - Fluid statistics simcomp, "fluid_stats", for more details see the [statistics guide](@ref statistics-guide) +- Computation of the spectral error indicator \ref simcomp_speri ## Controling execution and file output Each simulation component is, by default, executed once per time step to perform @@ -63,12 +64,25 @@ vorticity fields will be added to the main `.fld` file. ### vorticity {#simcomp_vorticity} Computes the vorticity field an stores in the field registry as `omega_x`, -`omega_y` and `omega_z`. +`omega_y` and `omega_z`. By default, appends the 3 vorticity fields to the field files as +scalars. To output in a different `fld` series, use the `"output_filename"` parameter. + + ~~~~~~~~~~~~~~~{.json} + { + "type": "vorticity" + } + ~~~~~~~~~~~~~~~ ### lambda2 {#simcomp_lambda2} Computes \f$ \lambda_2 \f$ for the velocity field and stores it in the normal output files as the first unused field. This means that \f$ \lambda_2 \f$ can be found in the temperature field in then fld files if running without a scalar -and s1 if neko is run with one scalar. +and s1 if neko is run with one scalar. To output in a different `fld` series, use the `"output_filename"` parameter. + + ~~~~~~~~~~~~~~~{.json} + { + "type": "lambda2" + } + ~~~~~~~~~~~~~~~ ### probes {#simcomp_probes} Probes selected solution fields at a list of points. This list of points can be @@ -198,15 +212,15 @@ field to derivate is controlled by the `field` keyword and the direction by the `direction` keyword. The simcomp will register the computed derivatives in the registry as `d[field]_d[direction]`, where the values in the brackets correspond to the choice of the user keywords. Supports writing the computed -fields to disk via the usual common keywords. +fields to disk via the usual common keywords. The resulting field will be +appended as a scalar to the field files. To output in a different `fld` series, +use the `"output_filename"` parameter. ~~~~~~~~~~~~~~~{.json} { "type": "derivative", "field": "u", - "direction", "y", - "output_control" : "simulation_time", - "output_value" : 1.0 + "direction": "y" } ~~~~~~~~~~~~~~~ @@ -247,3 +261,18 @@ writing the computed fields to disk via the usual common keywords. "output_control" : "never" } ~~~~~~~~~~~~~~~ + +### Spectral error indicator {#simcomp_speri} + +Computes the spectral error indicator as developed by Mavriplis (1989) (https://doi.org/10.1007/978-3-663-13975-1_34). +This is an a posteriori error measure, based on the local properties of +the spectral solution. This method formally only gives an indication of the error. + +The spectral error indicator is computed for the 3 velocity fields, resulting +in 3 additional fields appended to the field files. + +~~~~~~~~~~~~~~~{.json} + { + "type": "spectral_error" + } + ~~~~~~~~~~~~~~~ diff --git a/src/.depends b/src/.depends index c4a7b063ae4..524f43e9608 100644 --- a/src/.depends +++ b/src/.depends @@ -27,7 +27,6 @@ sem/map_2d.o : sem/map_2d.f90 field/field.o io/fld_file_data.o comm/mpi_types.o sem/dofmap.o : sem/dofmap.f90 mesh/hex.o mesh/quad.o mesh/element.o math/math.o device/device.o math/tensor.o math/fast3d.o common/utils.o config/num_types.o adt/tuple.o sem/space.o mesh/mesh.o config/neko_config.o sem/coef.o : sem/coef.f90 device/device.o math/mxm_wrapper.o sem/bcknd/device/device_coef.o math/bcknd/device/device_math.o mesh/mesh.o math/math.o sem/space.o sem/dofmap.o config/num_types.o config/neko_config.o gs/gs_ops.o gs/gather_scatter.o sem/cpr.o : sem/cpr.f90 sem/dofmap.o common/log.o math/mxm_wrapper.o math/tensor.o sem/coef.o mesh/mesh.o math/math.o sem/space.o field/field.o config/num_types.o config/neko_config.o gs/gather_scatter.o -sem/spectral_error_indicator.o : sem/spectral_error_indicator.f90 common/utils.o comm/comm.o device/device.o config/neko_config.o gs/gather_scatter.o math/bcknd/device/device_math.o math/tensor.o io/file.o math/math.o field/field_list.o sem/coef.o field/field.o config/num_types.o common/time_interpolator.o : common/time_interpolator.f90 math/fast3d.o common/utils.o math/math.o math/bcknd/device/device_math.o config/neko_config.o field/field.o config/num_types.o sem/interpolation.o : sem/interpolation.f90 sem/space.o math/bcknd/cpu/tensor_cpu.o math/tensor.o math/fast3d.o device/device.o config/num_types.o config/neko_config.o sem/point_interpolator.o : sem/point_interpolator.f90 config/neko_config.o math/bcknd/device/device_math.o device/device.o common/utils.o math/fast3d.o math/math.o mesh/point.o config/num_types.o sem/space.o math/tensor.o @@ -277,11 +276,12 @@ scalar/source_scalar.o : scalar/source_scalar.f90 math/bcknd/device/device_math. simulation_components/simulation_component.o : simulation_components/simulation_component.f90 common/json_utils.o common/time_based_controller.o case.o config/num_types.o simulation_components/simcomp_executor.o : simulation_components/simcomp_executor.f90 common/log.o common/utils.o case.o common/json_utils.o simulation_components/simulation_component.o config/num_types.o simulation_components/vorticity.o : simulation_components/vorticity.f90 simulation_components/field_writer.o common/json_utils.o io/fld_file_output.o case.o math/operators.o field/field.o field/field_registry.o simulation_components/simulation_component.o config/num_types.o +simulation_components/spectral_error.o : simulation_components/spectral_error.f90 field/field_registry.o case.o common/json_utils.o simulation_components/simulation_component.o simulation_components/field_writer.o common/utils.o comm/comm.o device/device.o common/log.o config/neko_config.o gs/gather_scatter.o math/bcknd/device/device_math.o math/tensor.o io/file.o math/math.o field/field_list.o sem/coef.o field/field.o config/num_types.o simulation_components/lambda2.o : simulation_components/lambda2.f90 device/device.o simulation_components/field_writer.o case.o math/operators.o field/field.o field/field_registry.o simulation_components/simulation_component.o config/num_types.o simulation_components/weak_grad.o : simulation_components/weak_grad.f90 simulation_components/field_writer.o common/json_utils.o io/fld_file_output.o case.o math/operators.o field/field.o field/field_registry.o simulation_components/simulation_component.o config/num_types.o simulation_components/derivative.o : simulation_components/derivative.f90 common/utils.o simulation_components/field_writer.o common/json_utils.o io/fld_file_output.o case.o math/operators.o field/field.o field/field_registry.o simulation_components/simulation_component.o config/num_types.o simulation_components/les_simcomp.o : simulation_components/les_simcomp.f90 simulation_components/field_writer.o common/json_utils.o les/les_model.o case.o simulation_components/simulation_component.o config/num_types.o -simulation_components/simulation_component_fctry.o : simulation_components/simulation_component_fctry.f90 simulation_components/derivative.o simulation_components/weak_grad.o simulation_components/field_writer.o common/utils.o simulation_components/les_simcomp.o simulation_components/probes.o simulation_components/lambda2.o simulation_components/fluid_stats_simcomp.o simulation_components/force_torque.o simulation_components/vorticity.o simulation_components/simulation_component.o +simulation_components/simulation_component_fctry.o : simulation_components/simulation_component_fctry.f90 simulation_components/spectral_error.o simulation_components/derivative.o simulation_components/weak_grad.o simulation_components/field_writer.o common/utils.o simulation_components/les_simcomp.o simulation_components/probes.o simulation_components/lambda2.o simulation_components/fluid_stats_simcomp.o simulation_components/force_torque.o simulation_components/vorticity.o simulation_components/simulation_component.o source_terms/source_term.o : source_terms/source_term.f90 field/field_list.o sem/coef.o config/num_types.o source_terms/coriolis_source_term.o : source_terms/coriolis_source_term.f90 source_terms/bcknd/cpu/coriolis_source_term_cpu.o common/utils.o config/neko_config.o sem/coef.o source_terms/source_term.o common/json_utils.o field/field_list.o config/num_types.o source_terms/bcknd/cpu/coriolis_source_term_cpu.o : source_terms/bcknd/cpu/coriolis_source_term_cpu.f90 field/field.o math/math.o field/field_list.o config/num_types.o @@ -316,5 +316,5 @@ wall_models/wall_model.o : wall_models/wall_model.f90 common/log.o comm/comm.o m wall_models/rough_log_law.o : wall_models/rough_log_law.f90 common/utils.o common/json_utils.o field/field_registry.o wall_models/wall_model.o config/neko_config.o sem/coef.o sem/dofmap.o config/num_types.o field/field.o wall_models/spalding.o : wall_models/spalding.f90 common/utils.o common/log.o common/json_utils.o field/field_registry.o wall_models/wall_model.o config/neko_config.o sem/coef.o sem/dofmap.o config/num_types.o field/field.o wall_models/wall_model_fctry.o : wall_models/wall_model_fctry.f90 common/json_utils.o common/utils.o wall_models/rough_log_law.o wall_models/spalding.o les/vreman.o wall_models/wall_model.o -neko.o : neko.f90 common/json_utils.o common/runtime_statistics.o bc/field_dirichlet_vector.o bc/field_dirichlet.o mesh/point_zone_registry.o mesh/point_zones/sphere_point_zone.o mesh/point_zones/box_point_zone.o mesh/point_zone.o sem/point_interpolator.o common/time_interpolator.o io/data_streamer.o simulation_components/simcomp_executor.o field/scratch_registry.o field/field_registry.o qoi/drag_torque.o common/system.o sem/spectral_error_indicator.o simulation_components/probes.o simulation_components/simulation_component.o math/tensor.o math/matrix.o math/vector.o scalar/scalar_user_source_term.o fluid/fluid_user_source_term.o field/field_list.o fluid/fluid_stats.o sem/cpr.o sem/map_2d.o sem/map_1d.o math/bcknd/device/device_math.o device/device.o common/jobctrl.o common/signal.o common/user_intf.o common/projection.o math/mathops.o math/operators.o simulation.o io/output.o io/output_controller.o case.o config/neko_config.o comm/parmetis.o math/ax.o bc/dirichlet.o bc/wall.o bc/bc.o sem/coef.o gs/gather_scatter.o comm/mpi_types.o field/field.o io/file.o common/global_interpolation.o math/mxm_wrapper.o io/format/map.o field/mesh_field.o mesh/point.o mesh/mesh.o adt/tuple.o adt/stack.o adt/uset.o adt/htable.o sem/space.o sem/dofmap.o sem/speclib.o math/math.o common/log.o common/utils.o comm/comm.o config/num_types.o +neko.o : neko.f90 common/json_utils.o common/runtime_statistics.o bc/field_dirichlet_vector.o bc/field_dirichlet.o mesh/point_zone_registry.o mesh/point_zones/sphere_point_zone.o mesh/point_zones/box_point_zone.o mesh/point_zone.o sem/point_interpolator.o common/time_interpolator.o io/data_streamer.o simulation_components/simcomp_executor.o field/scratch_registry.o field/field_registry.o qoi/drag_torque.o common/system.o simulation_components/spectral_error.o simulation_components/probes.o simulation_components/simulation_component.o math/tensor.o math/matrix.o math/vector.o scalar/scalar_user_source_term.o fluid/fluid_user_source_term.o field/field_list.o fluid/fluid_stats.o sem/cpr.o sem/map_2d.o sem/map_1d.o math/bcknd/device/device_math.o device/device.o common/jobctrl.o common/signal.o common/user_intf.o common/projection.o math/mathops.o math/operators.o simulation.o io/output.o io/output_controller.o case.o config/neko_config.o comm/parmetis.o math/ax.o bc/dirichlet.o bc/wall.o bc/bc.o sem/coef.o gs/gather_scatter.o comm/mpi_types.o field/field.o io/file.o common/global_interpolation.o math/mxm_wrapper.o io/format/map.o field/mesh_field.o mesh/point.o mesh/mesh.o adt/tuple.o adt/stack.o adt/uset.o adt/htable.o sem/space.o sem/dofmap.o sem/speclib.o math/math.o common/log.o common/utils.o comm/comm.o config/num_types.o driver.o : driver.f90 neko.o diff --git a/src/Makefile.am b/src/Makefile.am index 0693ba408a3..12b13b7bfd9 100644 --- a/src/Makefile.am +++ b/src/Makefile.am @@ -30,7 +30,6 @@ neko_fortran_SOURCES = \ sem/dofmap.f90\ sem/coef.f90\ sem/cpr.f90\ - sem/spectral_error_indicator.f90\ common/time_interpolator.f90\ sem/interpolation.f90\ sem/point_interpolator.f90\ @@ -280,6 +279,7 @@ neko_fortran_SOURCES = \ simulation_components/simulation_component.f90\ simulation_components/simcomp_executor.f90\ simulation_components/vorticity.f90\ + simulation_components/spectral_error.f90\ simulation_components/lambda2.f90\ simulation_components/weak_grad.f90\ simulation_components/derivative.f90\ diff --git a/src/neko.f90 b/src/neko.f90 index 7f790b30c4a..eab565561e0 100644 --- a/src/neko.f90 +++ b/src/neko.f90 @@ -102,7 +102,7 @@ module neko use simulation_component, only : simulation_component_t, & simulation_component_wrapper_t use probes, only : probes_t - use spectral_error_indicator + use spectral_error use system, only : system_cpu_name, system_cpuid use drag_torque, only : drag_torque_zone, drag_torque_facet, drag_torque_pt use field_registry, only : neko_field_registry diff --git a/src/simulation_components/simulation_component_fctry.f90 b/src/simulation_components/simulation_component_fctry.f90 index 5deef0cad4f..80e8508eada 100644 --- a/src/simulation_components/simulation_component_fctry.f90 +++ b/src/simulation_components/simulation_component_fctry.f90 @@ -43,16 +43,18 @@ use field_writer, only : field_writer_t use weak_grad, only : weak_grad_t use derivative, only : derivative_t + use spectral_error, only: spectral_error_t ! List of all possible types created by the factory routine - character(len=20) :: SIMCOMPS_KNOWN_TYPES(7) = [character(len=20) :: & + character(len=20) :: SIMCOMPS_KNOWN_TYPES(8) = [character(len=20) :: & "vorticity", & "lambda2", & "probes", & "les_model", & "field_writer", & "fluid_stats", & - "force_torque"] + "force_torque", & + "spectral_error"] contains @@ -92,6 +94,8 @@ module subroutine simulation_component_factory(object, json, case) allocate(force_torque_t::object) else if (trim(type_name) .eq. "fluid_stats") then allocate(fluid_stats_simcomp_t::object) + else if (trim(type_name) .eq. "spectral_error") then + allocate(spectral_error_t::object) else type_string = concat_string_array(SIMCOMPS_KNOWN_TYPES, & NEW_LINE('A') // "- ", .true.) diff --git a/src/sem/spectral_error_indicator.f90 b/src/simulation_components/spectral_error.f90 similarity index 80% rename from src/sem/spectral_error_indicator.f90 rename to src/simulation_components/spectral_error.f90 index 24805fae6f5..c57b7d98fba 100644 --- a/src/sem/spectral_error_indicator.f90 +++ b/src/simulation_components/spectral_error.f90 @@ -30,8 +30,8 @@ ! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE ! POSSIBILITY OF SUCH DAMAGE. ! -!> Implements type spectral_error_indicator_t. -module spectral_error_indicator +!> Implements type spectral_error_t. +module spectral_error use num_types, only: rp use field, only: field_t use coefs, only: coef_t @@ -42,9 +42,17 @@ module spectral_error_indicator use device_math, only: device_copy use gather_scatter use neko_config - use device, only: DEVICE_TO_HOST, device_memcpy + use logger, only: neko_log + use device, only: DEVICE_TO_HOST, HOST_TO_DEVICE, device_memcpy use comm, only: pe_rank - use utils, only: NEKO_FNAME_LEN + use utils, only: NEKO_FNAME_LEN, neko_error + use field_writer, only: field_writer_t + use simulation_component, only: simulation_component_t + use json_module, only: json_file + use json_utils, only: json_get, json_get_or_default + use case, only: case_t + use field_registry, only: neko_field_registry + use, intrinsic :: iso_c_binding implicit none private @@ -54,15 +62,15 @@ module spectral_error_indicator !! This is a posteriori error measure, based on the local properties of !! the spectral solution, which was developed by Mavriplis. This method !! formally only gives an indication of the error. - type, public :: spectral_error_indicator_t + type, public, extends(simulation_component_t) :: spectral_error_t !> Pointers to main fields type(field_t), pointer :: u => null() type(field_t), pointer :: v => null() type(field_t), pointer :: w => null() !> Transformed fields - type(field_t) :: u_hat - type(field_t) :: v_hat - type(field_t) :: w_hat + type(field_t), pointer :: u_hat => null() + type(field_t), pointer :: v_hat => null() + type(field_t), pointer :: w_hat => null() !> Working field - Consider making this a simple array type(field_t) :: wk !> Configuration of spectral error calculation @@ -86,46 +94,62 @@ module spectral_error_indicator type(field_list_t) :: speri_l !> File to write type(file_t) :: mf_speri + !> Field writer controller for the output + type(field_writer_t) :: writer contains - !> Initialize object. - procedure, pass(this) :: init => spec_err_ind_init - !> Destructor - procedure, pass(this) :: free => spec_err_ind_free - !> Calculate the indicator - procedure, pass(this) :: get_indicators => spec_err_ind_get - !> Calculate the indicator - procedure, pass(this) :: write_as_field => spec_err_ind_write + !> Constructor. + procedure, pass(this) :: init => spectral_error_init + !> Destructor. + procedure, pass(this) :: free => spectral_error_free + !> Compute the indicator (called according to the simcomp controller). + procedure, pass(this) :: compute_ => spectral_error_compute + !> Calculate the indicator. + procedure, pass(this) :: get_indicators => spectral_error_get_indicators - end type spectral_error_indicator_t + end type spectral_error_t contains - !> Constructor - !! @param u u velocity field - !! @param v v velocity field - !! @param w w velocity field - !! @param coef type with all geometrical variables - subroutine spec_err_ind_init(this, u,v,w,coef) - class(spectral_error_indicator_t), intent(inout) :: this - type(field_t), intent(in), target :: u - type(field_t), intent(in), target :: v - type(field_t), intent(in), target :: w + !> Constructor. + subroutine spectral_error_init(this, json, case) + class(spectral_error_t), intent(inout) :: this + type(json_file), intent(inout) :: json + class(case_t), intent(inout), target :: case + + character(len=20) :: fields(3) + + !> Add keyword "fields" to the json so that the field writer + ! picks it up. Will also add those fields to the registry. + fields(1) = "u_hat" + fields(2) = "v_hat" + fields(3) = "w_hat" + call json%add("fields", fields) + + call this%init_base(json, case) + call this%writer%init(json, case) + + call spectral_error_init_from_attributes(this, case%fluid%c_Xh) + + end subroutine spectral_error_init + + !> Actual constructor. + !! @param coef type with all geometrical variables. + subroutine spectral_error_init_from_attributes(this, coef) + class(spectral_error_t), intent(inout) :: this type(coef_t), intent(in) :: coef integer :: il, jl, aa character(len=NEKO_FNAME_LEN) :: fname_speri - !> call destructior - call this%free() + this%u => neko_field_registry%get_field("u") + this%v => neko_field_registry%get_field("v") + this%w => neko_field_registry%get_field("w") + this%u_hat => neko_field_registry%get_field("u_hat") + this%v_hat => neko_field_registry%get_field("v_hat") + this%w_hat => neko_field_registry%get_field("w_hat") - !> Assign the pointers - this%u => u - this%v => v - this%w => w !> Initialize fields and copy data from proper one - this%u_hat = u - this%v_hat = v - this%w_hat = w - this%wk = u + this%wk = this%u + !> Allocate arrays (Consider moving some to coef) allocate(this%eind_u(coef%msh%nelv)) allocate(this%eind_v(coef%msh%nelv)) @@ -134,7 +158,7 @@ subroutine spec_err_ind_init(this, u,v,w,coef) allocate(this%sig_v(coef%msh%nelv)) allocate(this%sig_w(coef%msh%nelv)) - !> The following code has been lifted from Adams implementation + !> The following code has been lifted from Adam's implementation associate(LX1 => coef%Xh%lx, LY1 => coef%Xh%ly, & LZ1 => coef%Xh%lz, & SERI_SMALL => this%SERI_SMALL, & @@ -147,29 +171,21 @@ subroutine spec_err_ind_init(this, u,v,w,coef) ) ! correctness check if (SERI_NP.gt.SERI_NP_MAX) then - if (pe_rank.eq.0) write(*,*) 'SETI_NP greater than SERI_NP_MAX' + call neko_log%message('SETI_NP greater than SERI_NP_MAX') endif il = SERI_NP+SERI_ELR jl = min(LX1,LY1) jl = min(jl,LZ1) if (il.gt.jl) then - if (pe_rank.eq.0) write(*,*) 'SERI_NP+SERI_ELR greater than L?1' + call neko_log%message('SERI_NP+SERI_ELR greater than L?1') endif end associate - !> Initialize the list that holds the fields to write - call list_init3(this%speri_l,this%u_hat,this%v_hat, & - this%w_hat) - ! Initialize the file - fname_speri = 'speri.fld' - this%mf_speri = file_t(fname_speri) - - end subroutine spec_err_ind_init - + end subroutine spectral_error_init_from_attributes !> Destructor - subroutine spec_err_ind_free(this) - class(spectral_error_indicator_t), intent(inout) :: this + subroutine spectral_error_free(this) + class(spectral_error_t), intent(inout) :: this if(allocated(this%eind_u)) then deallocate(this%eind_u) @@ -195,27 +211,65 @@ subroutine spec_err_ind_free(this) deallocate(this%sig_w) end if - call this%u_hat%free() - call this%v_hat%free() - call this%w_hat%free() call this%wk%free() nullify(this%u) nullify(this%v) nullify(this%w) + nullify(this%u_hat) + nullify(this%v_hat) + nullify(this%w_hat) + + call this%writer%free() + call this%free_base() + + end subroutine spectral_error_free + + !> Compute the spectral error indicator. + subroutine spectral_error_compute(this, t, tstep) + class(spectral_error_t), intent(inout) :: this + real(kind=rp), intent(in) :: t + integer, intent(in) :: tstep - !> finalize data related to writing - call list_final3(this%speri_l) - call file_free(this%mf_speri) + integer :: e, i, lx, ly, lz, nelv, n - end subroutine spec_err_ind_free - !> Transform a field u > u_hat into physical or spectral space - !! the result of the transformation is in u_hat - !! @param u_hat transformed field - !! @param u field to transform - !! @param wk working field - !! @param coef type coef for mesh parameters + call this%get_indicators(this%case%fluid%c_Xh) + + lx = this%u_hat%Xh%lx + ly = this%u_hat%Xh%ly + lz = this%u_hat%Xh%lz + nelv = this%u_hat%msh%nelv + + !> Copy the element indicator into all points of the field + do e = 1,nelv + do i = 1,lx*ly*lx + this%u_hat%x(i,1,1,e) = this%eind_u(e) + this%v_hat%x(i,1,1,e) = this%eind_v(e) + this%w_hat%x(i,1,1,e) = this%eind_w(e) + end do + end do + + ! We need this copy to GPU since the field writer does an internal copy + ! GPU -> CPU before writing the field files. This overwrites the *_hat + ! host arrays that contain the values. + if (NEKO_BCKND_DEVICE .eq. 1) then + call device_memcpy(this%u_hat%x, this%u_hat%x_d, lx*ly*lz*nelv, & + HOST_TO_DEVICE, sync = .true.) + call device_memcpy(this%v_hat%x, this%v_hat%x_d, lx*ly*lz*nelv, & + HOST_TO_DEVICE, sync = .true.) + call device_memcpy(this%w_hat%x, this%w_hat%x_d, lx*ly*lz*nelv, & + HOST_TO_DEVICE, sync = .true.) + end if + + end subroutine spectral_error_compute + + !> Transform a field u to u_hat into physical or spectral space + !! the result of the transformation is in u_hat. + !! @param u_hat Transformed field (output). + !! @param u Field to transform (input). + !! @param wk Working field. + !! @param coef Type coef for mesh parameters. !! @param space String that indicates which space to transform, "spec" or "phys". subroutine transform_to_spec_or_phys(u_hat, u, wk, coef, space) type(field_t), intent(inout) :: u_hat @@ -261,8 +315,8 @@ end subroutine transform_to_spec_or_phys !> Transform and get the spectral error indicators !! @param coef type coef for mesh parameters and space - subroutine spec_err_ind_get(this,coef) - class(spectral_error_indicator_t), intent(inout) :: this + subroutine spectral_error_get_indicators(this,coef) + class(spectral_error_t), intent(inout) :: this type(coef_t), intent(inout) :: coef integer :: i @@ -271,49 +325,34 @@ subroutine spec_err_ind_get(this,coef) call transform_to_spec_or_phys(this%v_hat, this%v, this%wk, coef, 'spec') call transform_to_spec_or_phys(this%w_hat, this%w, this%wk, coef, 'spec') - ! Get the spectral error indicator - call calculate_indicators(this, coef, this%eind_u, this%sig_u, coef%msh%nelv, & - coef%Xh%lx, coef%Xh%ly, coef%Xh%lz, & - this%u_hat%x) - call calculate_indicators(this, coef, this%eind_v, this%sig_v, coef%msh%nelv, & - coef%Xh%lx, coef%Xh%ly, coef%Xh%lz, & - this%v_hat%x) - call calculate_indicators(this, coef, this%eind_w, this%sig_w, coef%msh%nelv, & - coef%Xh%lx, coef%Xh%ly, coef%Xh%lz, & - this%w_hat%x) - - end subroutine spec_err_ind_get + call calculate_indicators(this, coef, this%eind_u, this%sig_u, & + coef%msh%nelv, coef%Xh%lx, coef%Xh%ly, coef%Xh%lz, & + this%u_hat%x) + call calculate_indicators(this, coef, this%eind_v, this%sig_v, & + coef%msh%nelv, coef%Xh%lx, coef%Xh%ly, coef%Xh%lz, & + this%v_hat%x) + call calculate_indicators(this, coef, this%eind_w, this%sig_w, & + coef%msh%nelv, coef%Xh%lx, coef%Xh%ly, coef%Xh%lz, & + this%w_hat%x) + + end subroutine spectral_error_get_indicators !> Write error indicators in a field file. !! @param t Current simulation time. - subroutine spec_err_ind_write(this, t) - class(spectral_error_indicator_t), intent(inout) :: this + subroutine spectral_error_write(this, t) + class(spectral_error_t), intent(inout) :: this real(kind=rp), intent(in) :: t integer i, e integer lx, ly, lz, nelv - lx = this%u_hat%Xh%lx - ly = this%u_hat%Xh%ly - lz = this%u_hat%Xh%lz - nelv = this%u_hat%msh%nelv - - !> Copy the element indicator into all points of the field - do e = 1,nelv - do i = 1,lx*ly*lx - this%u_hat%x(i,1,1,e) = this%eind_u(e) - this%v_hat%x(i,1,1,e) = this%eind_v(e) - this%w_hat%x(i,1,1,e) = this%eind_w(e) - end do - end do - !> Write the file !! Remember that the list is already ponting to the fields !! that were just modified. call this%mf_speri%write(this%speri_l,t) - end subroutine spec_err_ind_write + end subroutine spectral_error_write !> Wrapper for old fortran 77 subroutines !! @param coef coef type @@ -325,7 +364,7 @@ end subroutine spec_err_ind_write !! @param LZ1 gll points in z !! @paran var variable to calculate indicator subroutine calculate_indicators(this, coef, eind, sig, lnelt, LX1, LY1, LZ1, var) - type(spectral_error_indicator_t), intent(inout) :: this + type(spectral_error_t), intent(inout) :: this type(coef_t), intent(in) :: coef integer :: lnelt integer :: LX1 @@ -360,7 +399,7 @@ end subroutine calculate_indicators !! @param LY1 gll points in y !! @param LZ1 gll points in z subroutine speri_var(this, est,sig,var,nell,xa,xb,LX1,LY1,LZ1) - type(spectral_error_indicator_t), intent(inout) :: this + type(spectral_error_t), intent(inout) :: this integer :: nell integer :: LX1 integer :: LY1 @@ -472,7 +511,6 @@ subroutine speri_var(this, est,sig,var,nell,xa,xb,LX1,LY1,LZ1) end subroutine speri_var - !> Extrapolate the Legendre spectrum from the last points !! @param estx spectral indicator !! @param sigx coefficient of the exponential fit @@ -485,7 +523,7 @@ end subroutine speri_var subroutine speri_extrap(this,estx,sigx,coef11,coef, & ix_st,ix_en,nyl,nzl) implicit none - type(spectral_error_indicator_t), intent(inout) :: this + type(spectral_error_t), intent(inout) :: this !> argument list integer :: ix_st,ix_en,nyl,nzl !> Legendre coefficients; last SERI_NP columns @@ -664,33 +702,4 @@ subroutine speri_extrap(this,estx,sigx,coef11,coef, & end subroutine speri_extrap - !> Helper function to initialize field list to write - !! @param list list to allocate - !! @param uu field to add to the list - !! @param vv field to add to the list - !! @param ww field to add to the list - subroutine list_init3(list, uu, vv, ww) - type(field_list_t), intent(inout) :: list - type(field_t), target:: uu - type(field_t), target:: vv - type(field_t), target:: ww - !> Initialize field lists - call list%init(3) - call list%assign_to_field(1, uu) - call list%assign_to_field(2, vv) - call list%assign_to_field(3, ww) - end subroutine list_init3 - - - !> Helper function to finalize field list to write - !! @param list list to deallocate - subroutine list_final3(list) - type(field_list_t), intent(inout) :: list - call list%free - end subroutine list_final3 - - -end module spectral_error_indicator - - - +end module spectral_error