diff --git a/LibTrixi.jl/src/LibTrixi.jl b/LibTrixi.jl/src/LibTrixi.jl index 57d79722..ab855489 100644 --- a/LibTrixi.jl/src/LibTrixi.jl +++ b/LibTrixi.jl/src/LibTrixi.jl @@ -1,9 +1,9 @@ module LibTrixi using OrdinaryDiffEq: OrdinaryDiffEq, step!, check_error, DiscreteCallback -using Trixi: Trixi, summary_callback, mesh_equations_solver_cache, nelements, - nelementsglobal, nvariables, nnodes, wrap_array, eachelement, cons2prim, - get_node_vars, eachnode +using Trixi: Trixi, summary_callback, mesh_equations_solver_cache, ndims, nelements, + nelementsglobal, ndofs, ndofsglobal, nvariables, nnodes, wrap_array, + eachelement, cons2prim, get_node_vars, eachnode using MPI: MPI, run_init_hooks, set_default_error_handler_return using Pkg @@ -28,15 +28,36 @@ export trixi_ndims, export trixi_nelements, trixi_nelements_cfptr, trixi_nelements_jl -export trixi_nelements_global, - trixi_nelements_global_cfptr, - trixi_nelements_global_jl +export trixi_nelementsglobal, + trixi_nelementsglobal_cfptr, + trixi_nelementsglobal_jl +export trixi_ndofs, + trixi_ndofs_cfptr, + trixi_ndofs_jl +export trixi_ndofsglobal, + trixi_ndofsglobal_cfptr, + trixi_ndofsglobal_jl +export trixi_ndofselement, + trixi_ndofselement_cfptr, + trixi_ndofselement_jl export trixi_nvariables, trixi_nvariables_cfptr, trixi_nvariables_jl -export trixi_load_cell_averages, - trixi_load_cell_averages_cfptr, - trixi_load_cell_averages_jl +export trixi_nnodes, + trixi_nnodes_cfptr, + trixi_nnodes_jl +export trixi_load_node_reference_coordinates, + trixi_load_node_reference_coordinates_cfptr, + trixi_load_node_reference_coordinates_jl +export trixi_load_node_weights, + trixi_load_node_weights_cfptr, + trixi_load_node_weights_jl +export trixi_load_primitive_vars, + trixi_load_primitive_vars_cfptr, + trixi_load_primitive_vars_jl +export trixi_load_element_averaged_primitive_vars, + trixi_load_element_averaged_primitive_vars_cfptr, + trixi_load_element_averaged_primitive_vars_jl export trixi_version_library, trixi_version_library_cfptr, trixi_version_library_jl diff --git a/LibTrixi.jl/src/api_c.jl b/LibTrixi.jl/src/api_c.jl index e21ab685..3609ea24 100644 --- a/LibTrixi.jl/src/api_c.jl +++ b/LibTrixi.jl/src/api_c.jl @@ -285,7 +285,7 @@ trixi_ndims_cfptr() = @cfunction(trixi_ndims, Cint, (Cint,)) """ trixi_nelements(simstate_handle::Cint)::Cint -Return number of local elements (cells). +Return number of elements local to the MPI rank. """ function trixi_nelements end @@ -298,18 +298,63 @@ trixi_nelements_cfptr() = @cfunction(trixi_nelements, Cint, (Cint,)) """ - trixi_nelements_global(simstate_handle::Cint)::Cint + trixi_nelementsglobal(simstate_handle::Cint)::Cint -Return number of global elements (cells). +Return global number of elements on all MPI ranks. """ -function trixi_nelements_global end +function trixi_nelementsglobal end -Base.@ccallable function trixi_nelements_global(simstate_handle::Cint)::Cint +Base.@ccallable function trixi_nelementsglobal(simstate_handle::Cint)::Cint simstate = load_simstate(simstate_handle) - return trixi_nelements_global_jl(simstate) + return trixi_nelementsglobal_jl(simstate) end -trixi_nelements_global_cfptr() = @cfunction(trixi_nelements_global, Cint, (Cint,)) +trixi_nelementsglobal_cfptr() = @cfunction(trixi_nelementsglobal, Cint, (Cint,)) + + +""" + trixi_ndofs(simstate_handle::Cint)::Cint + +Return number of degrees of freedom (all quadrature points on all elements of current rank). +""" +function trixi_ndofs end + +Base.@ccallable function trixi_ndofs(simstate_handle::Cint)::Cint + simstate = load_simstate(simstate_handle) + return trixi_ndofs_jl(simstate) +end + +trixi_ndofs_cfptr() = @cfunction(trixi_ndofs, Cint, (Cint,)) + + +""" + trixi_ndofsglobal(simstate_handle::Cint)::Cint + +Return global number of degrees of freedom (all quadrature points on all elements on all ranks). +""" +function trixi_ndofsglobal end + +Base.@ccallable function trixi_ndofsglobal(simstate_handle::Cint)::Cint + simstate = load_simstate(simstate_handle) + return trixi_ndofsglobal_jl(simstate) +end + +trixi_ndofsglobal_cfptr() = @cfunction(trixi_ndofsglobal, Cint, (Cint,)) + + +""" + trixi_ndofselement(simstate_handle::Cint)::Cint + +Return number of degrees of freedom per element. +""" +function trixi_ndofselement end + +Base.@ccallable function trixi_ndofselement(simstate_handle::Cint)::Cint + simstate = load_simstate(simstate_handle) + return trixi_ndofselement_jl(simstate) +end + +trixi_ndofselement_cfptr() = @cfunction(trixi_ndofselement, Cint, (Cint,)) """ @@ -328,32 +373,122 @@ trixi_nvariables_cfptr() = @cfunction(trixi_nvariables, Cint, (Cint,)) """ - trixi_load_cell_averages(data::Ptr{Cdouble}, simstate_handle::Cint)::Cvoid + trixi_nnodes(simstate_handle::Cint)::Cint + +Return number of quadrature nodes per dimension. +""" +function trixi_nnodes end + +Base.@ccallable function trixi_nnodes(simstate_handle::Cint)::Cint + simstate = load_simstate(simstate_handle) + return trixi_nnodes_jl(simstate) +end + +trixi_nnodes_cfptr() = @cfunction(trixi_nnodes, Cint, (Cint,)) + + +""" + trixi_load_node_reference_coordinates(simstate_handle::Cint, data::Ptr{Cdouble})::Cvoid + +Get reference coordinates of 1D quadrature nodes. +""" +function trixi_load_node_reference_coordinates end + +Base.@ccallable function trixi_load_node_reference_coordinates(simstate_handle::Cint, + data::Ptr{Cdouble})::Cvoid + simstate = load_simstate(simstate_handle) + + # convert C to Julia array + size = trixi_nnodes_jl(simstate) + data_jl = unsafe_wrap(Array, data, size) + + trixi_load_node_reference_coordinates_jl(simstate, data_jl) + return nothing +end + +trixi_load_node_reference_coordinates_cfptr() = + @cfunction(trixi_load_node_reference_coordinates, Cvoid, (Cint, Ptr{Cdouble})) + + +""" + trixi_load_node_weights(simstate_handle::Cint, data::Ptr{Cdouble})::Cvoid + +Get weights of 1D quadrature nodes. +""" +function trixi_load_node_weights end + +Base.@ccallable function trixi_load_node_weights(simstate_handle::Cint, + data::Ptr{Cdouble})::Cvoid + simstate = load_simstate(simstate_handle) + + # convert C to Julia array + size = trixi_nnodes_jl(simstate) + data_jl = unsafe_wrap(Array, data, size) + + return trixi_load_node_weights_jl(simstate, data_jl) +end + +trixi_load_node_weights_cfptr() = + @cfunction(trixi_load_node_weights, Cvoid, (Cint, Ptr{Cdouble})) + + +""" + trixi_load_primitive_vars(simstate_handle::Cint, variable_id::Cint, + data::Ptr{Cdouble})::Cvoid + +Load primitive variable. + +The values for the primitive variable at position `variable_id` at every degree of freedom +are stored in the given array `data`. + +The given array has to be of correct size (ndofs) and memory has to be allocated beforehand. +""" +function trixi_load_primitive_vars end + +Base.@ccallable function trixi_load_primitive_vars(simstate_handle::Cint, variable_id::Cint, + data::Ptr{Cdouble})::Cvoid + simstate = load_simstate(simstate_handle) + + # convert C to Julia array + size = trixi_ndofs_jl(simstate) + data_jl = unsafe_wrap(Array, data, size) + + trixi_load_primitive_vars_jl(simstate, variable_id, data_jl) + return nothing +end + +trixi_load_primitive_vars_cfptr() = + @cfunction(trixi_load_primitive_vars, Cvoid, (Cint, Cint, Ptr{Cdouble})) + + +""" + trixi_load_element_averaged_primitive_vars(simstate_handle::Cint, variable_id::Cint, + data::Ptr{Cdouble})::Cvoid -Return cell averaged solution state. +Load element averages for primitive variable. -Cell averaged values for each cell and each primitive variable are stored in a contiguous -array, where cell values for the first variable appear first and values for the other -variables subsequently (structure-of-arrays layout). +Element averaged values for the primitive variable at position `variable_id` for each +element are stored in the given array `data`. -The given array has to be of correct size and memory has to be allocated beforehand. +The given array has to be of correct size (nelements) and memory has to be allocated +beforehand. """ -function trixi_load_cell_averages end +function trixi_load_element_averaged_primitive_vars end -Base.@ccallable function trixi_load_cell_averages(data::Ptr{Cdouble}, - simstate_handle::Cint)::Cvoid +Base.@ccallable function trixi_load_element_averaged_primitive_vars(simstate_handle::Cint, + variable_id::Cint, data::Ptr{Cdouble})::Cvoid simstate = load_simstate(simstate_handle) # convert C to Julia array - size = trixi_nvariables_jl(simstate) * trixi_nelements_jl(simstate) + size = trixi_nelements_jl(simstate) data_jl = unsafe_wrap(Array, data, size) - trixi_load_cell_averages_jl(data_jl, simstate) + trixi_load_element_averaged_primitive_vars_jl(simstate, variable_id, data_jl) return nothing end -trixi_load_cell_averages_cfptr() = - @cfunction(trixi_load_cell_averages, Cvoid, (Ptr{Cdouble}, Cint,)) +trixi_load_element_averaged_primitive_vars_cfptr() = + @cfunction(trixi_load_element_averaged_primitive_vars, Cvoid, (Cint, Cint, Ptr{Cdouble})) """ diff --git a/LibTrixi.jl/src/api_jl.jl b/LibTrixi.jl/src/api_jl.jl index 4eb6c70b..649420a5 100644 --- a/LibTrixi.jl/src/api_jl.jl +++ b/LibTrixi.jl/src/api_jl.jl @@ -86,22 +86,85 @@ function trixi_nelements_jl(simstate) end -function trixi_nelements_global_jl(simstate) +function trixi_nelementsglobal_jl(simstate) mesh, _, solver, cache = mesh_equations_solver_cache(simstate.semi) return nelementsglobal(mesh, solver, cache) end +function trixi_ndofs_jl(simstate) + mesh, _, solver, cache = mesh_equations_solver_cache(simstate.semi) + return ndofs(mesh, solver, cache) +end + + +function trixi_ndofsglobal_jl(simstate) + mesh, _, solver, cache = mesh_equations_solver_cache(simstate.semi) + return ndofsglobal(mesh, solver, cache) +end + + +function trixi_ndofselement_jl(simstate) + mesh, _, solver, _ = mesh_equations_solver_cache(simstate.semi) + return nnodes(solver)^ndims(mesh) +end + + function trixi_nvariables_jl(simstate) _, equations, _, _ = mesh_equations_solver_cache(simstate.semi) return nvariables(equations) end -function trixi_load_cell_averages_jl(data, simstate) +function trixi_nnodes_jl(simstate) + _, _, solver, _ = mesh_equations_solver_cache(simstate.semi) + return nnodes(solver) +end + + +function trixi_load_node_reference_coordinates_jl(simstate, data) + _, _, solver, _ = mesh_equations_solver_cache(simstate.semi) + for i in eachnode(solver) + data[i] = solver.basis.nodes[i] + end +end + + +function trixi_load_node_weights_jl(simstate, data) + _, _, solver, _ = mesh_equations_solver_cache(simstate.semi) + for i in eachnode(solver) + data[i] = solver.basis.weights[i] + end +end + + +function trixi_load_primitive_vars_jl(simstate, variable_id, data) + mesh, equations, solver, cache = mesh_equations_solver_cache(simstate.semi) + n_nodes_per_dim = nnodes(solver) + n_dims = ndims(mesh) + n_nodes = n_nodes_per_dim^n_dims + + u_ode = simstate.integrator.u + u = wrap_array(u_ode, mesh, equations, solver, cache) + + # all permutations of nodes indices for arbitrary dimension + node_cis = CartesianIndices(ntuple(i -> n_nodes_per_dim, n_dims)) + node_lis = LinearIndices(node_cis) + + for element in eachelement(solver, cache) + for node_ci in node_cis + node_vars = get_node_vars(u, equations, solver, node_ci, element) + node_index = (element-1) * n_nodes + node_lis[node_ci] + data[node_index] = cons2prim(node_vars, equations)[variable_id] + end + end + + return nothing +end + + +function trixi_load_element_averaged_primitive_vars_jl(simstate, variable_id, data) mesh, equations, solver, cache = mesh_equations_solver_cache(simstate.semi) - n_elements = nelements(solver, cache) - n_variables = nvariables(equations) n_nodes = nnodes(solver) n_dims = ndims(mesh) @@ -111,15 +174,13 @@ function trixi_load_cell_averages_jl(data, simstate) # all permutations of nodes indices for arbitrary dimension node_cis = CartesianIndices(ntuple(i -> n_nodes, n_dims)) - # temporary storage for mean value on current element for all variables - u_mean = get_node_vars(u, equations, solver, node_cis[1], 1) - for element in eachelement(solver, cache) # compute mean value using nodal dg values and quadrature - u_mean = zero(u_mean) + u_mean = zero(eltype(u)) for node_ci in node_cis - u_node_prim = cons2prim(get_node_vars(u, equations, solver, node_ci, element), equations) + u_node_prim = cons2prim(get_node_vars(u, equations, solver, node_ci, element), + equations)[variable_id] weight = 1. for node_index in Tuple(node_ci) weight *= solver.basis.weights[node_index] @@ -130,11 +191,8 @@ function trixi_load_cell_averages_jl(data, simstate) # normalize to unit element u_mean = u_mean / 2^n_dims - # copy to provided array - # all element values for first variable, then for second, ... - for ivar = 0:n_variables-1 - data[element + ivar * n_elements] = u_mean[ivar+1] - end + # write to provided array + data[element] = u_mean end return nothing diff --git a/LibTrixi.jl/test/test_interface.jl b/LibTrixi.jl/test/test_interface.jl index fd9f70a0..4cfec875 100644 --- a/LibTrixi.jl/test/test_interface.jl +++ b/LibTrixi.jl/test/test_interface.jl @@ -77,20 +77,59 @@ end nelements_jl = trixi_nelements_jl(simstate_jl) @test nelements_c == nelements_jl - nelements_global_c = trixi_nelements_global(handle) - nelements_global_jl = trixi_nelements_global_jl(simstate_jl) - @test nelements_global_c == nelements_global_jl + nelementsglobal_c = trixi_nelementsglobal(handle) + nelementsglobal_jl = trixi_nelementsglobal_jl(simstate_jl) + @test nelementsglobal_c == nelementsglobal_jl + + # compare number of dofs + ndofs_c = trixi_ndofs(handle) + ndofs_jl = trixi_ndofs_jl(simstate_jl) + @test ndofs_c == ndofs_jl + + ndofsglobal_c = trixi_ndofsglobal(handle) + ndofsglobal_jl = trixi_ndofsglobal_jl(simstate_jl) + @test ndofsglobal_c == ndofsglobal_jl + + ndofselement_c = trixi_ndofselement(handle) + ndofselement_jl = trixi_ndofselement_jl(simstate_jl) + @test ndofselement_c == ndofselement_jl # compare number of variables nvariables_c = trixi_nvariables(handle) nvariables_jl = trixi_nvariables_jl(simstate_jl) @test nvariables_c == nvariables_jl - # compare cell averaged values - data_c = zeros(nvariables_c * nelements_c) - trixi_load_cell_averages(pointer(data_c), handle) - data_jl = zeros(nvariables_jl * nelements_jl) - trixi_load_cell_averages_jl(data_jl, simstate_jl) + # compare number of quadrature nodes + nnodes_c = trixi_nnodes(handle) + nnodes_jl = trixi_nnodes_jl(simstate_jl) + @test nnodes_c == nnodes_jl + + # compare coordinates of quadrature nodes + data_c = zeros(nnodes_c) + trixi_load_node_reference_coordinates(handle, pointer(data_c)) + data_jl = zeros(nnodes_jl) + trixi_load_node_reference_coordinates_jl(simstate_jl, data_jl) + @test data_c == data_jl + + # compare weights of quadrature nodes + data_c = zeros(nnodes_c) + trixi_load_node_weights(handle, pointer(data_c)) + data_jl = zeros(nnodes_jl) + trixi_load_node_weights_jl(simstate_jl, data_jl) + @test data_c == data_jl + + # compare element averaged values + data_c = zeros(nelements_c) + trixi_load_element_averaged_primitive_vars(handle, Int32(1), pointer(data_c)) + data_jl = zeros(nelements_jl) + trixi_load_element_averaged_primitive_vars_jl(simstate_jl, 1, data_jl) + @test data_c == data_jl + + # compare primitive variable values on all dofs + data_c = zeros(ndofs_c) + trixi_load_primitive_vars(handle, Int32(1), pointer(data_c)) + data_jl = zeros(ndofs_jl) + trixi_load_primitive_vars_jl(simstate_jl, 1, data_jl) @test data_c == data_jl end diff --git a/README.md b/README.md index ac496087..744ba918 100644 --- a/README.md +++ b/README.md @@ -2,7 +2,7 @@ [![Docs-stable](https://img.shields.io/badge/docs-stable-blue.svg)](https://trixi-framework.github.io/libtrixi/stable) [![Docs-dev](https://img.shields.io/badge/docs-dev-blue.svg)](https://trixi-framework.github.io/libtrixi/dev) -[![Build Status](https://github.com/trixi-framework/libtrixi/workflows/CI/badge.svg)](https://github.com/trixi-framework/libtrixi/actions?query=workflow%3ACI) +[![Build Status](https://github.com/trixi-framework/libtrixi/actions/workflows/ci.yml/badge.svg)](https://github.com/trixi-framework/libtrixi/actions?query=workflow%3ACI) [![Coveralls](https://coveralls.io/repos/github/trixi-framework/libtrixi/badge.svg)](https://coveralls.io/github/trixi-framework/libtrixi) [![Codecov](https://codecov.io/gh/trixi-framework/libtrixi/branch/main/graph/badge.svg)](https://codecov.io/gh/trixi-framework/libtrixi) [![License: MIT](https://img.shields.io/badge/License-MIT-success.svg)](https://opensource.org/licenses/MIT) diff --git a/examples/trixi_controller_data.c b/examples/trixi_controller_data.c index 5e6d0f1f..35e88cdd 100644 --- a/examples/trixi_controller_data.c +++ b/examples/trixi_controller_data.c @@ -46,16 +46,16 @@ int main ( int argc, char *argv[] ) { printf("\n*** Trixi controller *** nelements %d\n", nelements); // Allocate memory - data = realloc( data, sizeof(double) * nelements * nvariables ); + data = realloc( data, sizeof(double) * nelements ); - // Get averaged cell values for each variable - trixi_load_cell_averages(data, handle); + // Get element averaged values for first variable + trixi_load_element_averaged_primitive_vars(handle, 1, data); } } // Print first variable for (int i = 0; i < nelements; ++i) { - printf("u[cell %3d] = %f\n", i, data[i]); + printf("u[element %3d] = %f\n", i, data[i]); } // Finalize Trixi simulation diff --git a/examples/trixi_controller_data.f90 b/examples/trixi_controller_data.f90 index cdb49ecf..d5d09667 100644 --- a/examples/trixi_controller_data.f90 +++ b/examples/trixi_controller_data.f90 @@ -63,16 +63,16 @@ program trixi_controller_data_f ! allocate memory if ( associated(data) ) deallocate(data) - allocate( data(nelements*nvariables) ) + allocate( data(nelements) ) - ! get averaged cell values for each variable - call trixi_load_cell_averages(data, handle) + ! get element averaged values for first variable + call trixi_load_element_averaged_primitive_vars(handle, 1, data) end if end do ! print first variable do i = 1,nelements - print "('u[cell ', i4, '] = ', e14.8)", i, data(i) + print "('u[element ', i4, '] = ', e14.8)", i, data(i) end do ! Finalize Trixi simulation diff --git a/examples/trixi_controller_mpi.c b/examples/trixi_controller_mpi.c index c778a048..aadf4213 100644 --- a/examples/trixi_controller_mpi.c +++ b/examples/trixi_controller_mpi.c @@ -12,7 +12,7 @@ void init_mpi_external ( int argc, char *argv[] ) { int flag_init; ret = MPI_Initialized(&flag_init); - printf("[EXT] MPI Initialized: return %d, initialized %d, MPI_COMM_WORLD %p\n", ret, flag_init, MPI_COMM_WORLD); + printf("[EXT] MPI Initialized: return %d, initialized %d, MPI_COMM_WORLD %p\n", ret, flag_init, (void *) MPI_COMM_WORLD); if ( flag_init == 0 ) { diff --git a/examples/trixi_controller_mpi.f90 b/examples/trixi_controller_mpi.f90 index 4c90f881..4e62ca39 100644 --- a/examples/trixi_controller_mpi.f90 +++ b/examples/trixi_controller_mpi.f90 @@ -9,9 +9,9 @@ subroutine init_mpi_external comm = MPI_COMM_WORLD call MPI_Initialized(flag_init, ierror) - write(*, '(a,i1,a,l,a,i0)') "[EXT] MPI Initialized: return ", ierror, & - ", initialized ", flag_init, & - ", MPI_COMM_WORLD ", comm + write(*, '(a,i1,a,l1,a,i0)') "[EXT] MPI Initialized: return ", ierror, & + ", initialized ", flag_init, & + ", MPI_COMM_WORLD ", comm if (.not.(flag_init)) then requested_threadlevel = MPI_THREAD_SERIALIZED diff --git a/src/api.c b/src/api.c index 0b9a0130..b6465c1a 100644 --- a/src/api.c +++ b/src/api.c @@ -15,8 +15,15 @@ enum { TRIXI_FTPR_NDIMS, TRIXI_FPTR_NELEMENTS, TRIXI_FPTR_NELEMENTS_GLOBAL, + TRIXI_FPTR_NDOFS, + TRIXI_FPTR_NDOFS_GLOBAL, + TRIXI_FPTR_NDOFS_ELEMENT, TRIXI_FTPR_NVARIABLES, - TRIXI_FTPR_LOAD_CELL_AVERAGES, + TRIXI_FPTR_NNODES, + TRIXI_FPTR_LOAD_NODE_REFERENCE_COORDINATES, + TRIXI_FPTR_LOAD_NODE_WEIGHTS, + TRIXI_FTPR_LOAD_PRIMITIVE_VARS, + TRIXI_FTPR_LOAD_ELEMENT_AVERAGED_PRIMITIVE_VARS, TRIXI_FTPR_VERSION_LIBRARY, TRIXI_FTPR_VERSION_LIBRARY_MAJOR, TRIXI_FTPR_VERSION_LIBRARY_MINOR, @@ -36,24 +43,31 @@ static void* trixi_function_pointers[TRIXI_NUM_FPTRS]; // List of function names to obtain C function pointers from Julia // OBS! If any name is longer than 250 characters, adjust buffer sizes in setup.c static const char* trixi_function_pointer_names[] = { - [TRIXI_FTPR_INITIALIZE_SIMULATION] = "trixi_initialize_simulation_cfptr", - [TRIXI_FTPR_CALCULATE_DT] = "trixi_calculate_dt_cfptr", - [TRIXI_FTPR_IS_FINISHED] = "trixi_is_finished_cfptr", - [TRIXI_FTPR_STEP] = "trixi_step_cfptr", - [TRIXI_FTPR_FINALIZE_SIMULATION] = "trixi_finalize_simulation_cfptr", - [TRIXI_FTPR_NDIMS] = "trixi_ndims_cfptr", - [TRIXI_FPTR_NELEMENTS] = "trixi_nelements_cfptr", - [TRIXI_FPTR_NELEMENTS_GLOBAL] = "trixi_nelements_global_cfptr", - [TRIXI_FTPR_NVARIABLES] = "trixi_nvariables_cfptr", - [TRIXI_FTPR_LOAD_CELL_AVERAGES] = "trixi_load_cell_averages_cfptr", - [TRIXI_FTPR_VERSION_LIBRARY] = "trixi_version_library_cfptr", - [TRIXI_FTPR_VERSION_LIBRARY_MAJOR] = "trixi_version_library_major_cfptr", - [TRIXI_FTPR_VERSION_LIBRARY_MINOR] = "trixi_version_library_minor_cfptr", - [TRIXI_FTPR_VERSION_LIBRARY_PATCH] = "trixi_version_library_patch_cfptr", - [TRIXI_FTPR_VERSION_JULIA] = "trixi_version_julia_cfptr", - [TRIXI_FTPR_VERSION_JULIA_EXTENDED] = "trixi_version_julia_extended_cfptr", - [TRIXI_FTPR_EVAL_JULIA] = "trixi_eval_julia_cfptr", - [TRIXI_FTPR_GET_T8CODE_FOREST] = "trixi_get_t8code_forest_cfptr" + [TRIXI_FTPR_INITIALIZE_SIMULATION] = "trixi_initialize_simulation_cfptr", + [TRIXI_FTPR_CALCULATE_DT] = "trixi_calculate_dt_cfptr", + [TRIXI_FTPR_IS_FINISHED] = "trixi_is_finished_cfptr", + [TRIXI_FTPR_STEP] = "trixi_step_cfptr", + [TRIXI_FTPR_FINALIZE_SIMULATION] = "trixi_finalize_simulation_cfptr", + [TRIXI_FTPR_NDIMS] = "trixi_ndims_cfptr", + [TRIXI_FPTR_NELEMENTS] = "trixi_nelements_cfptr", + [TRIXI_FPTR_NELEMENTS_GLOBAL] = "trixi_nelementsglobal_cfptr", + [TRIXI_FPTR_NDOFS] = "trixi_ndofs_cfptr", + [TRIXI_FPTR_NDOFS_GLOBAL] = "trixi_ndofsglobal_cfptr", + [TRIXI_FPTR_NDOFS_ELEMENT] = "trixi_ndofselement_cfptr", + [TRIXI_FTPR_NVARIABLES] = "trixi_nvariables_cfptr", + [TRIXI_FPTR_NNODES] = "trixi_nnodes_cfptr", + [TRIXI_FPTR_LOAD_NODE_REFERENCE_COORDINATES] = "trixi_load_node_reference_coordinates_cfptr", + [TRIXI_FPTR_LOAD_NODE_WEIGHTS] = "trixi_load_node_weights_cfptr", + [TRIXI_FTPR_LOAD_PRIMITIVE_VARS] = "trixi_load_primitive_vars_cfptr", + [TRIXI_FTPR_LOAD_ELEMENT_AVERAGED_PRIMITIVE_VARS] = "trixi_load_element_averaged_primitive_vars_cfptr", + [TRIXI_FTPR_VERSION_LIBRARY] = "trixi_version_library_cfptr", + [TRIXI_FTPR_VERSION_LIBRARY_MAJOR] = "trixi_version_library_major_cfptr", + [TRIXI_FTPR_VERSION_LIBRARY_MINOR] = "trixi_version_library_minor_cfptr", + [TRIXI_FTPR_VERSION_LIBRARY_PATCH] = "trixi_version_library_patch_cfptr", + [TRIXI_FTPR_VERSION_JULIA] = "trixi_version_julia_cfptr", + [TRIXI_FTPR_VERSION_JULIA_EXTENDED] = "trixi_version_julia_extended_cfptr", + [TRIXI_FTPR_EVAL_JULIA] = "trixi_eval_julia_cfptr", + [TRIXI_FTPR_GET_T8CODE_FOREST] = "trixi_get_t8code_forest_cfptr" }; // Track initialization/finalization status to prevent unhelpful errors @@ -449,13 +463,13 @@ int trixi_ndims(int handle) { /** * @anchor trixi_nelements_api_c * - * @brief Return number of local elements (cells). + * @brief Return number of local elements. * * These usually differ from the global count when doing parallel computations. * * @param[in] handle simulation handle * - * @see trixi_nelements_global_api_c + * @see trixi_nelementsglobal_api_c */ int trixi_nelements(int handle) { @@ -468,9 +482,9 @@ int trixi_nelements(int handle) { /** - * @anchor trixi_nelements_global_api_c + * @anchor trixi_nelementsglobal_api_c * - * @brief Return number of global elements (cells). + * @brief Return global number of elements. * * These usually differ from the local count when doing parallel computations. * @@ -478,13 +492,72 @@ int trixi_nelements(int handle) { * * @see trixi_nelements_api_c */ -int trixi_nelements_global(int handle) { +int trixi_nelementsglobal(int handle) { // Get function pointer - int (*nelements_global)(int) = trixi_function_pointers[TRIXI_FPTR_NELEMENTS_GLOBAL]; + int (*nelementsglobal)(int) = trixi_function_pointers[TRIXI_FPTR_NELEMENTS_GLOBAL]; // Call function - return nelements_global(handle); + return nelementsglobal(handle); +} + + +/** + * @anchor trixi_ndofs_api_c + * + * @brief Return number of local degrees of freedom. + * + * These usually differ from the global count when doing parallel computations. + * + * @param[in] handle simulation handle + * + * @see trixi_ndofsglobal_api_c + */ +int trixi_ndofs(int handle) { + + // Get function pointer + int (*ndofs)(int) = trixi_function_pointers[TRIXI_FPTR_NDOFS]; + + // Call function + return ndofs(handle); +} + + +/** + * @anchor trixi_ndofsglobal_api_c + * + * @brief Return global number of degrees of freedom. + * + * These usually differ from the local count when doing parallel computations. + * + * @param[in] handle simulation handle + * + * @see trixi_ndofs_api_c + */ +int trixi_ndofsglobal(int handle) { + + // Get function pointer + int (*ndofsglobal)(int) = trixi_function_pointers[TRIXI_FPTR_NDOFS_GLOBAL]; + + // Call function + return ndofsglobal(handle); +} + + +/** + * @anchor trixi_ndofselement_api_c + * + * @brief Return number of degrees of freedom per element. + * + * @param[in] handle simulation handle + */ +int trixi_ndofselement(int handle) { + + // Get function pointer + int (*ndofselement)(int) = trixi_function_pointers[TRIXI_FPTR_NDOFS_ELEMENT]; + + // Call function + return ndofselement(handle); } @@ -505,31 +578,118 @@ int trixi_nvariables(int handle) { } -/** - * @anchor trixi_load_cell_averages_api_c +/** + * @anchor trixi_nnodes_api_c + * + * @brief Return number of quadrature nodes per dimension. + * + * @param[in] handle simulation handle + */ +int trixi_nnodes(int handle) { + + // Get function pointer + int (*nnodes)(int) = trixi_function_pointers[TRIXI_FPTR_NNODES]; + + // Call function + return nnodes(handle); +} + + +/** + * @anchor trixi_load_node_reference_coordinates_api_c * - * @brief Return cell averaged values + * @brief Get reference coordinates of 1D quadrature nodes. * - * Cell averaged values for each cell and each primitive variable are stored in a - * contiguous array, where cell values for the first variable appear first and values for - * the other variables subsequently (structure-of-arrays layout). + * The reference coordinates in [-1,1] of the quadrature nodes in the current DG scheme are + * stored in the provided array `node_coords`. The given array has to be of correct size, + * i.e. `nnodes`, and memory has to be allocated beforehand. * - * The given array has to be of correct size and memory has to be allocated beforehand. + * @param[in] handle simulation handle + * @param[out] node_coords node reference coordinates + */ +void trixi_load_node_reference_coordinates(int handle, double* node_coords) { + + // Get function pointer + void (*load_node_reference_coordinates)(int, double *) = trixi_function_pointers[TRIXI_FPTR_LOAD_NODE_REFERENCE_COORDINATES]; + + // Call function + return load_node_reference_coordinates(handle, node_coords); +} + + +/** + * @anchor trixi_load_node_weights_api_c * - * @param[in] handle simulation handle - * @param[out] data cell averaged values for all cells and all primitive variables + * @brief Get weights of 1D quadrature nodes. + * + * The weights of the quadrature nodes in the current DG scheme are stored in the provided + * array `node_weights`. The given array has to be of correct size, i.e. `nnodes`, and + * memory has to be allocated beforehand. + * + * @param[in] handle simulation handle + * @param[out] node_weights node weights + */ +void trixi_load_node_weights(int handle, double* node_weights) { + + // Get function pointer + void (*load_node_weights)(int, double *) = trixi_function_pointers[TRIXI_FPTR_LOAD_NODE_WEIGHTS]; + + // Call function + return load_node_weights(handle, node_weights); +} + + +/** + * @anchor trixi_load_primitive_vars_api_c + * + * @brief Load primitive variable + * + * The values for the primitive variable at position `variable_id` at every degree of + * freedom are stored in the given array `data`. + * + * The given array has to be of correct size (ndofs) and memory has to be allocated + * beforehand. + * + * @param[in] handle simulation handle + * @param[in] variable_id index of variable + * @param[out] data values for all degrees of freedom */ -void trixi_load_cell_averages(double * data, int handle) { +void trixi_load_primitive_vars(int handle, int variable_id, double * data) { // Get function pointer - void (*load_cell_averages)(double *, int) = - trixi_function_pointers[TRIXI_FTPR_LOAD_CELL_AVERAGES]; + void (*load_primitive_vars)(int, int, double *) = + trixi_function_pointers[TRIXI_FTPR_LOAD_PRIMITIVE_VARS]; // Call function - load_cell_averages(data, handle); + load_primitive_vars(handle, variable_id, data); } +/** + * @anchor trixi_load_element_averaged_primitive_vars_api_c + * + * @brief Load element averages for primitive variable + * + * Element averaged values for the primitive variable at position `variable_id` for each + * element are stored in the given array `data`. + * + * The given array has to be of correct size (nelements) and memory has to be allocated + * beforehand. + * + * @param[in] handle simulation handle + * @param[in] variable_id index of variable + * @param[out] data element averaged values for all elements + */ +void trixi_load_element_averaged_primitive_vars(int handle, int variable_id, double * data) { + + // Get function pointer + void (*load_element_averaged_primitive_vars)(int, int, double *) = + trixi_function_pointers[TRIXI_FTPR_LOAD_ELEMENT_AVERAGED_PRIMITIVE_VARS]; + + // Call function + load_element_averaged_primitive_vars(handle, variable_id, data); +} + /******************************************************************************************/ /* T8code */ diff --git a/src/api.f90 b/src/api.f90 index 4ccd0929..6526bd0c 100644 --- a/src/api.f90 +++ b/src/api.f90 @@ -248,7 +248,7 @@ integer(c_int) function trixi_ndims(handle) bind(c) !> !! @fn LibTrixi::trixi_nelements::trixi_nelements(handle) !! - !! @brief Return number of local elements (cells) + !! @brief Return number of local elements !! !! @param[in] handle simulation handle !! @@ -259,14 +259,53 @@ integer(c_int) function trixi_nelements(handle) bind(c) end function !> - !! @fn LibTrixi::trixi_nelements_global::trixi_nelements_global(handle) + !! @fn LibTrixi::trixi_nelementsglobal::trixi_nelementsglobal(handle) !! - !! @brief Return number of global elements (cells) + !! @brief Return global number of elements !! !! @param[in] handle simulation handle !! - !! @see @ref trixi_nelements_global_api_c "trixi_nelements_global (C API)" - integer(c_int) function trixi_nelements_global(handle) bind(c) + !! @see @ref trixi_nelementsglobal_api_c "trixi_nelementsglobal (C API)" + integer(c_int) function trixi_nelementsglobal(handle) bind(c) + use, intrinsic :: iso_c_binding, only: c_int + integer(c_int), value, intent(in) :: handle + end function + + !> + !! @fn LibTrixi::trixi_ndofs::trixi_ndofs(handle) + !! + !! @brief Return number of local degrees of freedom + !! + !! @param[in] handle simulation handle + !! + !! @see @ref trixi_ndofs_api_c "trixi_ndofs (C API)" + integer(c_int) function trixi_ndofs(handle) bind(c) + use, intrinsic :: iso_c_binding, only: c_int + integer(c_int), value, intent(in) :: handle + end function + + !> + !! @fn LibTrixi::trixi_ndofsglobal::trixi_ndofsglobal(handle) + !! + !! @brief Return global number of degrees of freedom + !! + !! @param[in] handle simulation handle + !! + !! @see @ref trixi_ndofsglobal_api_c "trixi_ndofsglobal (C API)" + integer(c_int) function trixi_ndofsglobal(handle) bind(c) + use, intrinsic :: iso_c_binding, only: c_int + integer(c_int), value, intent(in) :: handle + end function + + !> + !! @fn LibTrixi::trixi_ndofselement::trixi_ndofselement(handle) + !! + !! @brief Return number of degrees of freedom per element. + !! + !! @param[in] handle simulation handle + !! + !! @see @ref trixi_ndofselement_api_c "trixi_ndofselement (C API)" + integer(c_int) function trixi_ndofselement(handle) bind(c) use, intrinsic :: iso_c_binding, only: c_int integer(c_int), value, intent(in) :: handle end function @@ -285,18 +324,88 @@ integer(c_int) function trixi_nvariables(handle) bind(c) end function !> - !! @fn LibTrixi::trixi_load_cell_averages::trixi_load_cell_averages(data, handle) + !! @fn LibTrixi::trixi_nnodes::trixi_nnodes(handle) !! - !! @brief Return cell averaged values + !! @brief Return number of quadrature nodes per dimension. !! !! @param[in] handle simulation handle - !! @param[out] data cell averaged values for all cells and all variables !! - !! @see @ref trixi_load_cell_averages_api_c "trixi_load_cell_averages (C API)" - subroutine trixi_load_cell_averages(data, handle) bind(c) + !! @see @ref trixi_nnodes_api_c "trixi_nnodes (C API)" + integer(c_int) function trixi_nnodes(handle) bind(c) + use, intrinsic :: iso_c_binding, only: c_int + integer(c_int), value, intent(in) :: handle + end function + + !> + !! @fn LibTrixi::trixi_load_node_reference_coordinates::trixi_load_node_reference_coordinates(handle, node_coords) + !! + !! @brief Get reference coordinates of 1D quadrature nodes. + !! + !! The reference coordinates in [-1,1] of the quadrature nodes in the current DG scheme are + !! stored in the provided array `node_coords`. The given array has to be of correct size, + !! i.e. `nnodes`, and memory has to be allocated beforehand. + !! + !! @param[in] handle simulation handle + !! @param[out] node_coords node reference coordinates + !! + !! @see @ref trixi_load_node_reference_coordinates_api_c "trixi_load_node_reference_coordinates (C API)" + subroutine trixi_load_node_reference_coordinates(handle, node_coords) bind(c) + use, intrinsic :: iso_c_binding, only: c_int, c_double + integer(c_int), value, intent(in) :: handle + real(c_double), dimension(*), intent(out) :: node_coords + end subroutine + + !> + !! @fn LibTrixi::trixi_load_node_weights::trixi_load_node_weights(handle, node_weights) + !! + !! @brief Get weights of 1D quadrature nodes. + !! + !! The weights of the quadrature nodes in the current DG scheme are stored in the provided + !! array `node_weights`. The given array has to be of correct size, i.e. `nnodes`, and + !! memory has to be allocated beforehand. + !! + !! @param[in] handle simulation handle + !! @param[out] node_weights node weights + !! + !! @see @ref trixi_load_node_weights_api_c "trixi_load_node_weights (C API)" + subroutine trixi_load_node_weights(handle, node_weights) bind(c) + use, intrinsic :: iso_c_binding, only: c_int, c_double + integer(c_int), value, intent(in) :: handle + real(c_double), dimension(*), intent(out) :: node_weights + end subroutine + + !> + !! @fn LibTrixi::trixi_load_primitive_vars::trixi_load_primitive_vars(handle, variable_id, data) + !! + !! @brief Load primitive variable + !! + !! @param[in] handle simulation handle + !! @param[in] variable_id index of variable + !! @param[out] data primitive variable values for all degrees of freedom + !! + !! @see @ref trixi_load_primitive_vars_api_c "trixi_load_primitive_vars (C API)" + subroutine trixi_load_primitive_vars(handle, variable_id, data) bind(c) use, intrinsic :: iso_c_binding, only: c_int, c_double + integer(c_int), value, intent(in) :: handle + integer(c_int), value, intent(in) :: variable_id real(c_double), dimension(*), intent(in) :: data + end subroutine + + !> + !! @fn LibTrixi::trixi_load_element_averaged_primitive_vars::trixi_load_element_averaged_primitive_vars(handle, variable_id, data) + !! + !! @brief Load element averages for primitive variable + !! + !! @param[in] handle simulation handle + !! @param[in] variable_id index of variable + !! @param[out] data averaged values for all elements + !! + !! @see @ref trixi_load_element_averaged_primitive_vars_api_c "trixi_load_element_averaged_primitive_vars (C API)" + subroutine trixi_load_element_averaged_primitive_vars(handle, variable_id, data) bind(c) + use, intrinsic :: iso_c_binding, only: c_int, c_double integer(c_int), value, intent(in) :: handle + integer(c_int), value, intent(in) :: variable_id + real(c_double), dimension(*), intent(in) :: data end subroutine diff --git a/src/trixi.h b/src/trixi.h index 9d8be9bd..93315148 100644 --- a/src/trixi.h +++ b/src/trixi.h @@ -27,10 +27,17 @@ void trixi_step(int handle); // Simulation data int trixi_ndims(int handle); int trixi_nelements(int handle); -int trixi_nelements_global(int handle); +int trixi_nelementsglobal(int handle); +int trixi_ndofs(int handle); +int trixi_ndofsglobal(int handle); +int trixi_ndofselement(int handle); int trixi_nvariables(int handle); +int trixi_nnodes(int handle); double trixi_calculate_dt(int handle); -void trixi_load_cell_averages(double * data, int handle); +void trixi_load_node_reference_coordinates(int handle, double* node_coords); +void trixi_load_node_weights(int handle, double* node_weights); +void trixi_load_primitive_vars(int handle, int variable_id, double * data); +void trixi_load_element_averaged_primitive_vars(int handle, int variable_id, double * data); // T8code #if !defined(T8_H) && !defined(T8_FOREST_GENERAL_H) diff --git a/test/c/simulation.cpp b/test/c/simulation.cpp index 0db32dc5..ff732d1e 100644 --- a/test/c/simulation.cpp +++ b/test/c/simulation.cpp @@ -59,63 +59,103 @@ TEST(CInterfaceTest, SimulationRun) { // Check number of elements int nelements = trixi_nelements(handle); - int nelements_global = trixi_nelements_global(handle); - EXPECT_EQ(nelements * nranks, nelements_global); + int nelementsglobal = trixi_nelementsglobal(handle); + EXPECT_EQ(nelements * nranks, nelementsglobal); + + // Check number of dofs + int ndofs = trixi_ndofs(handle); + int ndofsglobal = trixi_ndofsglobal(handle); + EXPECT_EQ(ndofs * nranks, ndofsglobal); + + int ndofselement = trixi_ndofselement(handle); + EXPECT_EQ(nelements * ndofselement, ndofs); + EXPECT_EQ(nelementsglobal * ndofselement, ndofsglobal); // Check number of variables int nvariables = trixi_nvariables(handle); EXPECT_EQ(nvariables, 4); - // Check cell averaged values - int size = nelements * nvariables; - std::vector cell_averages(size); - trixi_load_cell_averages(cell_averages.data(), handle); + // Check number of quadrature nodes + int nnodes = trixi_nnodes(handle); + EXPECT_EQ(nnodes, 5); + + // Check quadrature, integrate f(x) = x^4 over [-1,1] + std::vector nodes(nnodes); + std::vector weights(nnodes); + trixi_load_node_reference_coordinates(handle, nodes.data()); + trixi_load_node_weights(handle, weights.data()); + double integral = 0.0; + for (int i = 0; i < nnodes; ++i) { + integral += weights[i] * nodes[i] * nodes[i] * nodes[i]* nodes[i]; + } + EXPECT_NEAR(integral, 0.4, 1e-17); + + // Check primitive variable values on all dofs + std::vector rho(ndofs); + std::vector energy(ndofs); + trixi_load_primitive_vars(handle, 1, rho.data()); + trixi_load_primitive_vars(handle, 4, energy.data()); + // check memory boarders + EXPECT_DOUBLE_EQ(rho[0], 1.0); + EXPECT_DOUBLE_EQ(rho[ndofs-1], 1.0); + EXPECT_DOUBLE_EQ(energy[0], 1.0e-5); + EXPECT_DOUBLE_EQ(energy[ndofs-1], 1.0e-5); + + // Check element averaged values + std::vector rho_averages(nelements); + std::vector v1_averages(nelements); + std::vector v2_averages(nelements); + std::vector e_averages(nelements); + trixi_load_element_averaged_primitive_vars(handle, 1, rho_averages.data()); + trixi_load_element_averaged_primitive_vars(handle, 2, v1_averages.data()); + trixi_load_element_averaged_primitive_vars(handle, 3, v2_averages.data()); + trixi_load_element_averaged_primitive_vars(handle, 4, e_averages.data()); if (nranks == 1) { // check memory boarders (densities at the beginning, energies at the end) - EXPECT_DOUBLE_EQ(cell_averages[0], 1.0); - EXPECT_DOUBLE_EQ(cell_averages[size-1], 1.0e-5); + EXPECT_DOUBLE_EQ(rho_averages[0], 1.0); + EXPECT_DOUBLE_EQ(e_averages[nelements-1], 1.0e-5); // check values somewhere near the center (expect symmetries) // densities - EXPECT_NEAR(cell_averages[93], 0.88263491354796, 1e-14); - EXPECT_NEAR(cell_averages[94], 0.88263491354796, 1e-14); - EXPECT_NEAR(cell_averages[161], 0.88263491354796, 1e-14); - EXPECT_NEAR(cell_averages[162], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[ 93], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[ 94], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[161], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[162], 0.88263491354796, 1e-14); // velocities - EXPECT_NEAR(cell_averages[ 93+ nelements], -0.14037267400591, 1e-14); - EXPECT_NEAR(cell_averages[ 94+2*nelements], -0.14037267400591, 1e-14); - EXPECT_NEAR(cell_averages[161+2*nelements], 0.14037267400591, 1e-14); - EXPECT_NEAR(cell_averages[162+ nelements], 0.14037267400591, 1e-14); + EXPECT_NEAR(v1_averages[ 93], -0.14037267400591, 1e-14); + EXPECT_NEAR(v2_averages[ 94], -0.14037267400591, 1e-14); + EXPECT_NEAR(v2_averages[161], 0.14037267400591, 1e-14); + EXPECT_NEAR(v1_averages[162], 0.14037267400591, 1e-14); } else if (nranks == 2) { if (rank == 0) { // check memory boarders (densities at the beginning, energies at the end) - EXPECT_DOUBLE_EQ(cell_averages[0], 1.0); - EXPECT_DOUBLE_EQ(cell_averages[size-1], 1.0e-5); + EXPECT_DOUBLE_EQ(rho_averages[0], 1.0); + EXPECT_DOUBLE_EQ(e_averages[nelements-1], 1.0e-5); // check values somewhere near the center (expect symmetries) // densities - EXPECT_NEAR(cell_averages[93], 0.88263491354796, 1e-14); - EXPECT_NEAR(cell_averages[94], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[93], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[94], 0.88263491354796, 1e-14); // velocities - EXPECT_NEAR(cell_averages[93+ nelements], -0.14037267400591, 1e-14); - EXPECT_NEAR(cell_averages[94+2*nelements], -0.14037267400591, 1e-14); + EXPECT_NEAR(v1_averages[93], -0.14037267400591, 1e-14); + EXPECT_NEAR(v2_averages[94], -0.14037267400591, 1e-14); } else { // check memory boarders (densities at the beginning, energies at the end) - EXPECT_DOUBLE_EQ(cell_averages[0], 1.0); - EXPECT_DOUBLE_EQ(cell_averages[size-1], 1.0e-5); + EXPECT_DOUBLE_EQ(rho_averages[0], 1.0); + EXPECT_DOUBLE_EQ(e_averages[nelements-1], 1.0e-5); // check values somewhere near the center (expect symmetries) // densities - EXPECT_NEAR(cell_averages[33], 0.88263491354796, 1e-14); - EXPECT_NEAR(cell_averages[34], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[33], 0.88263491354796, 1e-14); + EXPECT_NEAR(rho_averages[34], 0.88263491354796, 1e-14); // velocities - EXPECT_NEAR(cell_averages[33+2*nelements], 0.14037267400591, 1e-14); - EXPECT_NEAR(cell_averages[34+ nelements], 0.14037267400591, 1e-14); + EXPECT_NEAR(v2_averages[33], 0.14037267400591, 1e-14); + EXPECT_NEAR(v1_averages[34], 0.14037267400591, 1e-14); } } else { FAIL() << "Test cannot be run with " << nranks << " ranks."; } - + // Finalize Trixi simulation trixi_finalize_simulation(handle); diff --git a/test/fortran/simulationRun_suite.f90 b/test/fortran/simulationRun_suite.f90 index 1f81eb6f..9dc5e57d 100644 --- a/test/fortran/simulationRun_suite.f90 +++ b/test/fortran/simulationRun_suite.f90 @@ -22,12 +22,13 @@ end subroutine collect_simulationRun_suite subroutine test_simulationRun(error) type(error_type), allocatable, intent(out) :: error - integer :: handle, ndims, nelements, nelements_global, nvariables, size + integer :: handle, ndims, nelements, nelementsglobal, nvariables, ndofsglobal, & + ndofselement, ndofs, size, nnodes, i logical :: finished_status ! dp as defined in test-drive integer, parameter :: dp = selected_real_kind(15) - real(dp) :: dt - real(dp), dimension(:), allocatable :: data + real(dp) :: dt, integral + real(dp), dimension(:), allocatable :: data, weights ! Initialize Trixi call trixi_initialize(julia_project_path) @@ -55,21 +56,58 @@ subroutine test_simulationRun(error) nelements = trixi_nelements(handle) call check(error, nelements, 256) - nelements_global = trixi_nelements_global(handle) - call check(error, nelements_global, 256) + nelementsglobal = trixi_nelementsglobal(handle) + call check(error, nelementsglobal, 256) + + ! Check number of dofs + ndofselement = trixi_ndofselement(handle) + call check(error, ndofselement, 25) + + ndofs = trixi_ndofs(handle) + call check(error, ndofs, nelements * ndofselement) + + ndofsglobal = trixi_ndofsglobal(handle) + call check(error, ndofsglobal, nelementsglobal * ndofselement) ! Check number of variables nvariables = trixi_nvariables(handle) call check(error, nvariables, 4) - ! Check cell averaged values - size = nelements*nvariables + ! Check number of quadrature nodes + nnodes = trixi_nnodes(handle) + call check(error, nnodes, 5) + + ! Check quadrature, integrate f(x) = x^4 over [-1,1] + size = nnodes allocate(data(size)) - call trixi_load_cell_averages(data, handle) - call check(error, data(1), 1.0_dp) - call check(error, data(929), 2.6605289164377273_dp) - call check(error, data(size), 1e-5_dp) - + allocate(weights(size)) + call trixi_load_node_reference_coordinates(handle, data) + call trixi_load_node_weights(handle, weights) + integral = 0.0_dp + do i = 1, size + integral = integral + weights(i) * data(i) * data(i) * data(i)* data(i) + end do + call check(error, integral, 0.4_dp) + deallocate(data) + + ! Check primitive variable values + size = ndofs + allocate(data(size)) + call trixi_load_primitive_vars(handle, 1, data) + call check(error, data(1), 1.0_dp) + call check(error, data(3200), 1.0_dp) + call check(error, data(size), 1.0_dp) + deallocate(data) + + ! Check element averaged values + size = nelements + allocate(data(size)) + call trixi_load_element_averaged_primitive_vars(handle, 1, data) + call check(error, data(1), 1.0_dp) + call check(error, data(94), 0.99833232379996562_dp) + call check(error, data(size), 1.0_dp) + deallocate(data) + ! Finalize Trixi simulation call trixi_finalize_simulation(handle)