Skip to content

LakeInitialize 中 lake/element 下标错用与入湖河段下标约定不一致 #10

Description

@DankerMu

问题概述

在检查当前 master3aec657)时,发现 LakeInitialize() 中有一处明确的 lake/element 下标错用,以及同一段 lake 初始化逻辑里更严重的 0/1-based 下标约定不一致。两者都发生在湖泊相关初始化路径中。

1. lake[i].area 累加时误用了 Ele[i]

位置:

  • src/ModelData/MD_Lake.cpp:62-65

当前逻辑:

for (int j = 0; j < NumEle ; j++){
    if(Ele[j].iLake == i + 1){ /* Lake element */
        lake[i].NumEleLake++;
        lake[i].area += Ele[i].area;
    }
}

这里外层 i 是 lake index,内层 j 才是 element index。判断条件使用的是 Ele[j].iLake == i + 1,但面积累加却用了 Ele[i].area

预期应为:

lake[i].area += Ele[j].area;

影响

这会把湖泊面积算成“湖内单元数量 × 第 i 个 element 的面积”,而不是湖内所有 element 的面积总和。

当前看起来 lake[i].area 没有参与主求解路径;湖泊水位变化实际使用的是 y2LakeArea[i],来自 bathymetry 的 lake[i].u_toparea。因此这更像是确定的元数据/未来输出风险,而不一定解释当前 solver 数值差异。但它是明确的 i/j 下标错用,应修复。

2. Riv::toLake 计数和填充存在 0/1-based 混用

位置:

  • src/ModelData/MD_Lake.cpp:50
  • src/ModelData/MD_Lake.cpp:79-83
  • src/ModelData/MD_Lake.cpp:123-126
  • 相关使用:src/ModelData/MD_RiverFlux.cpp:17-24

Riv[i].toLake 初始化为 0-based:

Riv[i].toLake = (-3 - Riv[i].down) - 1;  /* 0 ~ inf */

计数 NumRivIn 时也按 0-based 使用:

if(Riv[j].toLake == i){
    lake[i].NumRivIn++;
}

solver 里 river-to-lake flux 也按 0-based 访问:

QLakeRivIn[Riv[i].toLake] += QrivDown[i];

但填充 lake[i].iRivIn 时用了 i + 1

if(Riv[j].toLake == i + 1){
    lake[i].iRivIn[ii++] = j;
}

影响

这会导致 lake[i].iRivIn 列表错位或漏填。更严重的是,数组分配大小来自 Riv[j].toLake == i 的计数,但写入来自 Riv[j].toLake == i + 1 的数量;如果不同 lake 的入湖河段数量不同,可能发生堆越界写。

当前主求解的 QLakeRivIn 直接使用 Riv[i].toLake,不依赖 lake[i].iRivIn,因此主要风险在初始化数据结构一致性和内存安全。但这仍然是真实 bug。

frLake/iRivOut 也有类似计数 == i、填充 == i + 1 的模式。当前 frLake 似乎只初始化为 -1,未看到有效赋值来源;修复时建议一并统一 index convention。

附带发现:Element::slope 分母疑似 typo,但当前无 solver 影响

位置:

  • src/classes/Element.cpp:124-126

当前逻辑:

slope[0] =  (z_surf - pz1)/ d1;
slope[1] =  (z_surf - pz2)/ d1;
slope[2] =  (z_surf - pz3)/ d1;

d1/d2/d3 分别对应三个 edge 的 centroid-to-edge 距离,因此 slope[1]slope[2] 理论上应分别除以 d2d3

不过当前检索显示 slope[] 字段没有参与主 solver 通量计算,实际通量使用的是 Dist2Edge / Dist2Nabor。所以这是低优先级 typo/清理项,不建议把它作为当前水文数值偏差来源。

建议修复

  1. lake[i].area += Ele[i].area 改为 lake[i].area += Ele[j].area
  2. 统一 toLake/frLakeNumRivIn/NumRivOutiRivIn/iRivOut 的下标约定。按当前代码注释和 solver 使用,toLake 应保持 0-based,则填充处应使用 Riv[j].toLake == i
  3. ii/jj 写入增加边界断言或一致性检查,防止计数和填充条件再次漂移。
  4. 可选:修复 Element::slope[1/2] 分母 typo,并在 PR 中明确说明当前无 solver 影响。

建议验证

  • 构造最小 lake case:至少两个 lake,且每个 lake 入湖 river 数量不同。
  • 验证 NumRivIniRivIn 内容一致。
  • 使用 ASan 跑初始化路径,确认无 heap-buffer-overflow。
  • 对含 lake 的短模拟做 smoke test,确认 solver 输出不崩溃。

Metadata

Metadata

Assignees

No one assigned

    Labels

    bugSomething isn't working

    Type

    No type

    Fields

    No fields configured for issues without a type.

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions