`SolutionArray.getState(loc)` returns the state of its `Solution` object rather than the state found at `loc`

[LiemDQ]

Problem description

In the C++ API, calling SolutionArray.getState(loc) returns the state of the associated Solution object rather than the state contained at loc if loc was the last-accessed location in the array.

Steps to reproduce

#include "cantera/core.h"
#include "cantera/base/SolutionArray.h"
#include <cassert>

int main() {
    using namespace Cantera;
    
    int loc = 0;
    shared_ptr<Solution> sol = newSolution("h2o2.yaml", "ohmech");
    shared_ptr<SolutionArray> array = SolutionArray::create(sol, 1);
    vector<double> initial_state(sol->thermo()->stateSize());
    
    sol->thermo()->setTemperature(100.0);
    sol->thermo()->saveState(initial_state);
    array->updateState(loc); //state @ 0 has temperature of 100.0 K
    
    vector<double> second_state(sol->thermo()->stateSize());
    sol->thermo()->setTemperature(234.0);
    sol->thermo()->saveState(second_state);

    vector<double> retrieved_state = array->getState(loc); //should still be the state @ 0

    assert(retrieved_state[0] == initial_state[0]);
    assert(retrieved_state[0] != second_state[0]);

    return 0;
}

Behavior

Expected: getState(loc) should always return the state stored at loc, regardless of subsequent modifications to the associated Solution object.

In this case, the first entry of the state (temperature) is equal to 234.0, but it should be equal to 100.0.

System information

• Cantera version: 3.2.0
• Python/MATLAB/other software versions: C++

Additional context

The problem seems to stem from the interaction of getState with setLoc. getState first uses setLoc to copy the state data to the associated Solution object (by setting restore = true), and then copies the state of the Solution object to the output vector. The problem is that setLoc has a check to see whether the "active" index m_loc has changed, and skips the copying step entirely if it hasn't.

cantera/src/base/SolutionArray.cpp

Lines 813 to 822 in ef09fb4

	} else if (static_cast<size_t>(m_active[loc_]) == m_loc) {
	return;
	} else if (loc_ >= m_size) {
	throw IndexError("SolutionArray::setLoc", "indices", loc_, m_size);
	}
	m_loc = static_cast<size_t>(m_active[loc_]);
	if (restore) {
	size_t nState = m_sol->thermo()->stateSize();
	m_sol->thermo()->restoreState(nState, m_data->data() + m_loc * m_stride);
	}

The simplest fix, in my opinion, would be to change getState to directly copy the state data to the output vector:

vector<double> SolutionArray::getState(int loc)
{
    setLoc(loc, false); //update m_loc but don't restore Solution state
    size_t nState = m_sol->thermo()->stateSize();
    double* data_ptr = m_data->data() + m_loc * m_stride;
    vector<double> out(data_ptr, data_ptr + nState); //copy data directly to vector using iterator constructor
    return out;
}

I would also consider removing the static_cast<size_t>(m_active[loc_]) == m_loc optimization from setLoc. It doesn't meaningfully improve performance if restore is false and when restore is true, the implicit state coupling can lead to subtle edge cases like this one.