Hacktoberfest 2026: los issues que los mantenedores marcaron para octubre, abiertos y aptos para principiantes. Explorar issues de Hacktoberfest

Possible upper/lower triangle storage inconsistency in rates() in multi-phase equilibrium solver

Abierto
  • ¿Algo no está bien?
#331 0 comentarios 0 reacciones 0 asignados Ver en GitHub

Nadie ha tomado este issue todavía.

Evaluación

Dificultad
3/5
Tiempo estimado
Medio día
Aptitud para principiantes
50/100
Tipo de issue
Error
Claridad
Bastante claro
Estado de actividad
Activo
Stack tecnológico
cpp
Área
backend

Línea de trabajo

The suspect code is in rates() in src/thermo/MultiPhaseEquilSolver.cpp: the loop assembling A(m,p) only fills the lower triangle, and the solve near line 1136 uses selfadjointView<Upper>(). Read the assembly loop and that solve to confirm whether the upper triangle is filled elsewhere before the solve. Done means either switching to selfadjointView<Lower>() or populating the upper triangle, with a maintainer confirming which one is intended and the solver results checked against a known equilibrium case.

Escrito por el modelo de indexación a partir del texto del issue.

Descripción

I believe there may be an inconsistency in the rates() routine between the way the matrix A is assembled and the triangular part subsequently used by Eigen's self-adjoint solver.

A is populated with:

for (int p = 0; p < npr; ++p) {
    for (int m = p; m < npr; ++m) {
        for (int j = p_sizes[m]; j < p_sizes[m+1]; ++j) {
            sum = 0.0;
            for (int i = 0; i < ncr; ++i)
                sum += H(j,i)*dlamy(i,p);
            A(m,p) += y[j]*sum;
        }
    }
}

Since m >= p, this explicitly computes the lower triangular part of A, i.e. A(m,p) for m >= p. This appears intentional, since A is expected to be symmetric.

However, the resulting system is then solved with
https://github.com/mutationpp/Mutationpp/blob/7607ef4d67f4de115705f78f467b1ea48bf02c26/src/thermo/MultiPhaseEquilSolver.cpp#L1136

As I understand Eigen's SelfAdjointView, selfadjointView<Upper>() treats the upper triangular portion of the matrix as the stored representation of the self-adjoint matrix and does not use the corresponding entries in the lower triangle. I believe this is therefore causing the off-diagonal elements to be discarded in the calculation of dx.

In other words, since the assembly produces

[ a  0  0
  b  c  0
  d  e  f ]

the intended symmetric matrix is presumably

[ a  b  d
  b  c  e
  d  e  f ]

However, since the upper-triangular off-diagonal elements remain zero, selfadjointView<Upper>() represents the matrix as

[ a  0  0
  0  c  0
  0  0  f ]

instead.

Should selfadjointView<Lower>() be used here instead? Alternatively, is there some reason the upper triangle is expected to be populated elsewhere before this solve?

Lenguaje dominante
C++
Estrellas
145
Forks
79
Merge medio
1 d 15 h
PR fusionados (30 d)
2

Preparar el entorno

Primeros pasos

  1. Lee el issue completo y luego la guía de contribución del proyecto.
  2. Comenta en el issue que vas a ocuparte — evita que dos personas hagan lo mismo.
  3. Haz un fork del repositorio y trabaja en una rama.
  4. Abre un pull request que haga referencia al número del issue.

Más de mutationpp/Mutationpp

Todos los issues de mutationpp/Mutationpp

Issues similares

Más issues de C++

Recibe los nuevos issues en tu correo

Un resumen breve de issues de GitHub para principiantes.