diff --git a/include/user_functions/ZalesakSphereMassFlowRateKernel.h b/include/user_functions/ZalesakSphereMassFlowRateKernel.h new file mode 100644 index 000000000..9df161e12 --- /dev/null +++ b/include/user_functions/ZalesakSphereMassFlowRateKernel.h @@ -0,0 +1,38 @@ +// Copyright 2017 National Technology & Engineering Solutions of Sandia, LLC +// (NTESS), National Renewable Energy Laboratory, University of Texas Austin, +// Northwest Research Associates. Under the terms of Contract DE-NA0003525 +// with NTESS, the U.S. Government retains certain rights in this software. +// +// This software is released under the BSD 3-clause license. See LICENSE file +// for more details. +// + +#ifndef ZALESAKSPHEREMASSFLOWRATEALG_H +#define ZALESAKSPHEREMASSFLOWRATEALG_H + +#include "AssembleEdgeSolverAlgorithm.h" + +namespace sierra { +namespace nalu { + +class ZalesakSphereMassFlowRateEdgeAlg : public AssembleEdgeSolverAlgorithm +{ +public: + // TODO: refactor to use FieldManager + ZalesakSphereMassFlowRateEdgeAlg( + Realm&, stk::mesh::Part*, EquationSystem*, const bool = false); + + virtual ~ZalesakSphereMassFlowRateEdgeAlg() = default; + + virtual void execute(); + +private: + unsigned coordinates_{stk::mesh::InvalidOrdinal}; + unsigned edgeAreaVec_{stk::mesh::InvalidOrdinal}; + unsigned massFlowRate_{stk::mesh::InvalidOrdinal}; +}; + +} // namespace nalu +} // namespace sierra + +#endif /* ZALESAKSPHEREMASSFLOWRATEALG_H */ diff --git a/include/user_functions/ZalesakSphereVOFAuxFunction.h b/include/user_functions/ZalesakSphereVOFAuxFunction.h new file mode 100644 index 000000000..8f9a6f05e --- /dev/null +++ b/include/user_functions/ZalesakSphereVOFAuxFunction.h @@ -0,0 +1,42 @@ +// Copyright 2017 National Technology & Engineering Solutions of Sandia, LLC +// (NTESS), National Renewable Energy Laboratory, University of Texas Austin, +// Northwest Research Associates. Under the terms of Contract DE-NA0003525 +// with NTESS, the U.S. Government retains certain rights in this software. +// +// This software is released under the BSD 3-clause license. See LICENSE file +// for more details. +// + +#ifndef ZalesakSphereVOFAuxFunction_h +#define ZalesakSphereVOFAuxFunction_h + +#include + +#include + +namespace sierra { +namespace nalu { + +class ZalesakSphereVOFAuxFunction : public AuxFunction +{ +public: + ZalesakSphereVOFAuxFunction(); + + virtual ~ZalesakSphereVOFAuxFunction() {} + + using AuxFunction::do_evaluate; + virtual void do_evaluate( + const double* coords, + const double time, + const unsigned spatialDimension, + const unsigned numPoints, + double* fieldPtr, + const unsigned fieldSize, + const unsigned beginPos, + const unsigned endPos) const; +}; + +} // namespace nalu +} // namespace sierra + +#endif diff --git a/src/VolumeOfFluidEquationSystem.C b/src/VolumeOfFluidEquationSystem.C index 0339dc531..8058180e8 100644 --- a/src/VolumeOfFluidEquationSystem.C +++ b/src/VolumeOfFluidEquationSystem.C @@ -29,6 +29,7 @@ // edge kernels #include "edge_kernels/VOFAdvectionEdgeAlg.h" #include "user_functions/ZalesakDiskMassFlowRateKernel.h" +#include "user_functions/ZalesakSphereMassFlowRateKernel.h" // node kernels #include "node_kernels/NodeKernelUtils.h" @@ -40,6 +41,7 @@ #include "ngp_algorithms/NodalGradBndryElemAlg.h" #include "stk_topology/topology.hpp" #include "user_functions/ZalesakDiskVOFAuxFunction.h" +#include "user_functions/ZalesakSphereVOFAuxFunction.h" #include "user_functions/DropletVOFAuxFunction.h" #include "ngp_utils/NgpFieldBLAS.h" #include "ngp_utils/NgpLoopUtils.h" @@ -491,6 +493,27 @@ VolumeOfFluidEquationSystem::register_initial_condition_fcn( new ZalesakDiskMassFlowRateEdgeAlg(realm_, part, this, useAvgMdot); realm_.initCondAlg_.push_back(VOFSetMassFlowRate); } + } else if (fcnName == "zalesak_sphere") { + theAuxFunc = new ZalesakSphereVOFAuxFunction(); + // Initialize mass flow rate until momentum connection implemented + { + const bool useAvgMdot = (realm_.solutionOptions_->turbulenceModel_ == + TurbulenceModel::SST_AMS) + ? true + : false; + ScalarFieldType* density_ = + realm_.meta_data().get_field( + stk::topology::NODE_RANK, "density"); + std::vector userSpec(1); + userSpec[0] = 1.0; + AuxFunction* constantAuxFunc = new ConstantAuxFunction(0, 1, userSpec); + AuxFunctionAlgorithm* constantAuxAlg = new AuxFunctionAlgorithm( + realm_, part, density_, constantAuxFunc, stk::topology::NODE_RANK); + realm_.initCondAlg_.push_back(constantAuxAlg); + auto VOFSetMassFlowRate = + new ZalesakSphereMassFlowRateEdgeAlg(realm_, part, this, useAvgMdot); + realm_.initCondAlg_.push_back(VOFSetMassFlowRate); + } } else if (fcnName == "droplet") { theAuxFunc = new DropletVOFAuxFunction(); } else { diff --git a/src/user_functions/CMakeLists.txt b/src/user_functions/CMakeLists.txt index 1f767794b..fa3efa46a 100644 --- a/src/user_functions/CMakeLists.txt +++ b/src/user_functions/CMakeLists.txt @@ -37,5 +37,7 @@ target_sources(nalu PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/WindEnergyTaylorVortexPressureAuxFunction.C ${CMAKE_CURRENT_SOURCE_DIR}/ZalesakDiskVOFAuxFunction.C ${CMAKE_CURRENT_SOURCE_DIR}/ZalesakDiskMassFlowRateKernel.C + ${CMAKE_CURRENT_SOURCE_DIR}/ZalesakSphereVOFAuxFunction.C + ${CMAKE_CURRENT_SOURCE_DIR}/ZalesakSphereMassFlowRateKernel.C ${CMAKE_CURRENT_SOURCE_DIR}/DropletVOFAuxFunction.C ) diff --git a/src/user_functions/ZalesakSphereMassFlowRateKernel.C b/src/user_functions/ZalesakSphereMassFlowRateKernel.C new file mode 100644 index 000000000..85c7c00d6 --- /dev/null +++ b/src/user_functions/ZalesakSphereMassFlowRateKernel.C @@ -0,0 +1,79 @@ +// Copyright 2017 National Technology & Engineering Solutions of Sandia, LLC +// (NTESS), National Renewable Energy Laboratory, University of Texas Austin, +// Northwest Research Associates. Under the terms of Contract DE-NA0003525 +// with NTESS, the U.S. Government retains certain rights in this software. +// +// This software is released under the BSD 3-clause license. See LICENSE file +// for more details. +// + +#include "user_functions/ZalesakSphereMassFlowRateKernel.h" +#include "EquationSystem.h" +#include "PecletFunction.h" +#include "SolutionOptions.h" +#include "utils/StkHelpers.h" +#include "edge_kernels/EdgeKernelUtils.h" +#include "stk_mesh/base/NgpField.hpp" +#include "stk_mesh/base/Types.hpp" + +namespace sierra { +namespace nalu { + +ZalesakSphereMassFlowRateEdgeAlg::ZalesakSphereMassFlowRateEdgeAlg( + Realm& realm, + stk::mesh::Part* part, + EquationSystem* eqSystem, + const bool useAverages) + : AssembleEdgeSolverAlgorithm(realm, part, eqSystem) +{ + const auto& meta = realm.meta_data(); + + coordinates_ = get_field_ordinal(meta, realm.get_coordinates_name()); + edgeAreaVec_ = + get_field_ordinal(meta, "edge_area_vector", stk::topology::EDGE_RANK); + massFlowRate_ = get_field_ordinal( + meta, (useAverages) ? "average_mass_flow_rate" : "mass_flow_rate", + stk::topology::EDGE_RANK); +} + +void +ZalesakSphereMassFlowRateEdgeAlg::execute() +{ + const double eps = 1.0e-16; + const int ndim = realm_.meta_data().spatial_dimension(); + + // Defaults from AMR-Wind + const double period = 6.0; + const double xrot = 0.5; + const double yrot = 0.5; + + // STK stk::mesh::NgpField instances for capture by lambda + const auto& fieldMgr = realm_.ngp_field_manager(); + const auto coordinates = fieldMgr.get_field(coordinates_); + const auto edgeAreaVec = fieldMgr.get_field(edgeAreaVec_); + const auto massFlowRate = fieldMgr.get_field(massFlowRate_); + + run_algorithm( + realm_.bulk_data(), + KOKKOS_LAMBDA( + ShmemDataType & smdata, const stk::mesh::FastMeshIndex& edge, + const stk::mesh::FastMeshIndex& nodeL, + const stk::mesh::FastMeshIndex& nodeR) { + // Scratch work array for edgeAreaVector + NALU_ALIGNED DblType av[NDimMax_]; + // Populate area vector work array + for (int d = 0; d < ndim; ++d) + av[d] = edgeAreaVec.get(edge, d); + + NALU_ALIGNED DblType edge_centroid[NDimMax_]; + for (int d = 0; d < ndim; ++d) + edge_centroid[d] = + 0.5 * coordinates.get(nodeL, d) + 0.5 * coordinates.get(nodeR, d); + + massFlowRate.get(edge, 0) = 2.0 * M_PI / period * + ((yrot - edge_centroid[1]) * av[0] + (edge_centroid[0] - xrot) * av[1]); + }); +} + +} // namespace nalu +} // namespace sierra diff --git a/src/user_functions/ZalesakSphereVOFAuxFunction.C b/src/user_functions/ZalesakSphereVOFAuxFunction.C new file mode 100644 index 000000000..5103f3724 --- /dev/null +++ b/src/user_functions/ZalesakSphereVOFAuxFunction.C @@ -0,0 +1,66 @@ +// Copyright 2017 National Technology & Engineering Solutions of Sandia, LLC +// (NTESS), National Renewable Energy Laboratory, University of Texas Austin, +// Northwest Research Associates. Under the terms of Contract DE-NA0003525 +// with NTESS, the U.S. Government retains certain rights in this software. +// +// This software is released under the BSD 3-clause license. See LICENSE file +// for more details. +// + +#include +#include + +// basic c++ +#include +#include +#include + +namespace sierra { +namespace nalu { + +ZalesakSphereVOFAuxFunction::ZalesakSphereVOFAuxFunction() : AuxFunction(0, 1) +{ + // does nothing +} + +void +ZalesakSphereVOFAuxFunction::do_evaluate( + const double* coords, + const double /*time*/, + const unsigned spatialDimension, + const unsigned numPoints, + double* fieldPtr, + const unsigned fieldSize, + const unsigned /*beginPos*/, + const unsigned /*endPos*/) const +{ + for (unsigned p = 0; p < numPoints; ++p) { + + const double x = coords[0]; + const double y = coords[1]; + const double z = coords[2]; + + // These are the default arguments from the corresponding case in AMR-Wind + const double xc = 0.5; + const double yc = 0.72; + const double zc = 0.24; + const double radius = 0.16; + const double depth = 0.2; + const double width = 0.04; + + fieldPtr[0] = 0.0; + // Put VOF in sphere + if ((x - xc) * (x - xc) + (y - yc) * (y - yc) + (z - zc) * (z - zc) < radius * radius) + fieldPtr[0] = 1.0; + + // Remove slot + if (x - xc > -0.5 * width && x - xc < 0.5 * width && y - yc > radius - depth) + fieldPtr[0] = 0.0; + + fieldPtr += fieldSize; + coords += spatialDimension; + } +} + +} // namespace nalu +} // namespace sierra