HydStructures.hpp

Go to the documentation of this file.

// SPDX-License-Identifier: Apache-2.0
//
// Copyright 2026 Caleb Buahin
//
// Licensed under the Apache License, Version 2.0 (the "License");
// you may not use this file except in compliance with the License.
// You may obtain a copy of the License at
//
//     http://www.apache.org/licenses/LICENSE-2.0
//
// Unless required by applicable law or agreed to in writing, software
// distributed under the License is distributed on an "AS IS" BASIS,
// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
// See the License for the specific language governing permissions and
// limitations under the License.

 
#ifndef OPENSWMM_HYD_STRUCTURES_HPP
#define OPENSWMM_HYD_STRUCTURES_HPP
 
#include "../data/LinkData.hpp"
#include "../data/NodeData.hpp"
#include "../data/TableData.hpp"
#include <vector>
 
namespace openswmm {
 
struct SimulationContext;
 
namespace hydstruct {
 
// ============================================================================
// Constants
// ============================================================================
 
constexpr double GRAVITY = 32.2;
constexpr double FUDGE   = 0.0001;
 
// ============================================================================
// Flap gates
// ============================================================================

bool flapGateBlocks(const SimulationContext& ctx, int j, int n1, int n2, double q);
 
// ============================================================================
// Per-type SoA groups (like XSectGroups but for structure types)
// ============================================================================
 
struct PumpGroup {
    int count = 0;
    std::vector<int>    link_idx;      
    std::vector<int>    curve_idx;     
    std::vector<int>    curve_type;    
    std::vector<double> speed;         
    std::vector<double> y_on;          
    std::vector<double> y_off;         
    void resize(int n);
};
 
struct OrificeGroup {
    int count = 0;
    std::vector<int>    link_idx;
    std::vector<int>    shape;         
    // No cached cOrif/cWeir/hCrit: legacy orifice_setSetting recomputes them
    // from the CURRENT setting every step, so computeOrificeFlowK derives
    // them inline per call (see the legacy-exact block there) — an init-time
    // snapshot at full opening could only ever go stale.
    // No cached flap-gate flag: the kernels read ctx.links.has_flap_gate (and
    // flapGateBlocks) directly, so a stale init-time snapshot cannot diverge.
    std::vector<double> surf_area;
    std::vector<double> length_eff;
    void resize(int n);
};
 
struct WeirGroup {
    int count = 0;
    std::vector<int>    link_idx;
    std::vector<int>    weir_type;     
    std::vector<double> c_disch1;      
    std::vector<double> c_disch2;      
    std::vector<double> end_con;       
    std::vector<double> slope;         
    std::vector<int>    cd_curve;      
    // No cached flap-gate flag — see the note on OrificeGroup.
    std::vector<double> surf_area;
    std::vector<double> length_eff;
    std::vector<uint8_t> can_surcharge;
    void resize(int n);
};
 
struct OutletGroup {
    int count = 0;
    std::vector<int>    link_idx;
    std::vector<int>    curve_idx;     
    std::vector<double> q_coeff;       
    std::vector<double> q_expon;       
    void resize(int n);
};
 
// ============================================================================
// Structure flow solver
// ============================================================================
 
class StructureSolver {
public:
    void init(SimulationContext& ctx);

    void computeAllFlows(SimulationContext& ctx, double dt,
                         double* node_new_surf_area = nullptr);

    void updatePumpTargetSettings(SimulationContext& ctx);

    const std::vector<int>& nonConduitIndices() const noexcept { return nc_indices_; }

    void computeNonConduitFlowOne(SimulationContext& ctx, double dt,
                                  double* node_new_surf_area, int link_idx);

    void scatterHeldSurfArea(SimulationContext& ctx,
                             double* node_new_surf_area, int link_idx);
 
private:
    PumpGroup    pumps_;
    OrificeGroup orifices_;
    WeirGroup    weirs_;
    OutletGroup  outlets_;
    std::vector<int> nc_indices_;   
    std::vector<int> nc_group_k_;   

    void computePumpFlows(SimulationContext& ctx, double dt,
                          const double* node_new_surf_area);
 
    // Per-element bodies of the batch loops (k = row in the type group).
    void computePumpFlowK(SimulationContext& ctx, double dt,
                          const double* node_new_surf_area, int k);
    void computeOrificeFlowK(SimulationContext& ctx,
                             double* node_new_surf_area, int k);
    void computeWeirFlowK(SimulationContext& ctx,
                          double* node_new_surf_area, int k);
    void computeOutletFlowK(SimulationContext& ctx, int k);
    void scatterOrificeSurfArea(SimulationContext& ctx,
                                double* node_new_surf_area,
                                std::size_t uk, std::size_t uj,
                                std::size_t un1, std::size_t un2);

    void computeOrificeFlows(SimulationContext& ctx,
                             double* node_new_surf_area = nullptr);

    void computeWeirFlows(SimulationContext& ctx,
                          double* node_new_surf_area = nullptr);

    void computeOutletFlows(SimulationContext& ctx);
};
 
} // namespace hydstruct
} // namespace openswmm
 
#endif // OPENSWMM_HYD_STRUCTURES_HPP