ChemicalState.hpp

// Reaktoro is a unified framework for modeling chemically reactive systems.
//
// Copyright © 2014-2024 Allan Leal
//
// This library is free software; you can redistribute it and/or
// modify it under the terms of the GNU Lesser General Public
// License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
//
// This library is distributed in the hope that it will be useful,
// but WITHOUT ANY WARRANTY; without even the implied warranty of
// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
// Lesser General Public License for more details.
//
// You should have received a copy of the GNU Lesser General Public License
// along with this library. If not, see <http://www.gnu.org/licenses/>.
 
#pragma once
 
// Reaktoro includes
#include <Reaktoro/Common/Matrix.hpp>
#include <Reaktoro/Common/Types.hpp>
#include <Reaktoro/Core/ChemicalProps.hpp>
#include <Reaktoro/Core/ChemicalSystem.hpp>
 
// Forward declarations (Optima)
namespace Optima { class State; }
 
namespace Reaktoro {
 
//=================================================================================================
//
// ChemicalState
//
//=================================================================================================
 
class ChemicalState
{
public:
    // Forward declarations
    class Equilibrium;
 
    explicit ChemicalState(ChemicalSystem const& system);
 
    ChemicalState(ChemicalState const& other);
 
    virtual ~ChemicalState();
 
    auto operator=(ChemicalState other) -> ChemicalState&;
 
    // --------------------------------------------------------------------------------------------
    // METHODS FOR SETTING/GETTING TEMPERATURE
    // --------------------------------------------------------------------------------------------
 
    auto setTemperature(real const& value) -> void;
 
    auto setTemperature(real value, Chars unit) -> void;
 
    auto temperature(real const& value) -> void;
 
    auto temperature(real value, Chars unit) -> void;
 
    auto temperature() const -> real;
 
    // --------------------------------------------------------------------------------------------
    // METHODS FOR SETTING/GETTING PRESSURE
    // --------------------------------------------------------------------------------------------
 
    auto setPressure(real const& value) -> void;
 
    auto setPressure(real value, Chars unit) -> void;
 
    auto pressure(real const& value) -> void;
 
    auto pressure(real value, Chars unit) -> void;
 
    auto pressure() const -> real;
 
    // --------------------------------------------------------------------------------------------
    // METHODS FOR SETTING THE AMOUNT OR MASS OF SPECIES
    // --------------------------------------------------------------------------------------------
 
    auto setSpeciesAmounts(real const& value) -> void;
 
    auto setSpeciesAmounts(ArrayXrConstRef const& n) -> void;
 
    auto setSpeciesAmounts(ArrayXdConstRef const& n) -> void;
 
    auto setSpeciesAmount(Index ispecies, real const& amount) -> void;
 
    auto setSpeciesAmount(StringOrIndex const& species, real amount, Chars unit) -> void;
 
    auto setSpeciesMass(StringOrIndex const& species, real mass, Chars unit) -> void;
 
    auto set(StringOrIndex const& species, real value, Chars unit) -> void;
 
    auto add(StringOrIndex const& species, real value, Chars unit) -> void;
 
    // --------------------------------------------------------------------------------------------
    // METHODS FOR GETTING THE AMOUNT OR MASS OF SPECIES, ELEMENTS, AND CHARGE
    // --------------------------------------------------------------------------------------------
 
    auto speciesAmounts() const -> ArrayXrConstRef;
 
    auto speciesAmountsInPhase(StringOrIndex const& phase) const -> ArrayXrConstRef;
 
    auto speciesAmount(StringOrIndex const& species) const -> real;
 
    auto speciesMass(StringOrIndex const& species) const -> real;
 
    auto componentAmounts() const -> ArrayXr;
 
    auto elementAmounts() const -> ArrayXr;
 
    auto charge() const -> real;
 
    // --------------------------------------------------------------------------------------------
    // METHODS TO SCALE THE AMOUNTS OF SPECIES IN THE SYSTEM OR PART OF IT
    // --------------------------------------------------------------------------------------------
 
    auto scaleSpeciesAmounts(real const& scalar) -> void;
 
    auto scaleSpeciesAmounts(real const& scalar, Indices const& indices) -> void;
 
    auto scaleSpeciesAmountsInPhase(StringOrIndex const& phase, real const& scalar) -> void;
 
    // --------------------------------------------------------------------------------------------
    // METHODS TO SCALE THE AMOUNT OF THE SYSTEM OR PART OF IT
    // --------------------------------------------------------------------------------------------
 
    auto scaleAmount(real amount, Chars unit) -> void;
 
    auto scalePhaseAmount(StringOrIndex const& phase, real amount, Chars unit) -> void;
 
    auto scaleFluidAmount(real amount, Chars unit) -> void;
 
    auto scaleSolidAmount(real amount, Chars unit) -> void;
 
    // --------------------------------------------------------------------------------------------
    // METHODS TO SCALE THE MASS OF THE SYSTEM OR PART OF IT
    // --------------------------------------------------------------------------------------------
 
    auto scaleMass(real mass, Chars unit) -> void;
 
    auto scalePhaseMass(StringOrIndex const& phase, real mass, Chars unit) -> void;
 
    auto scaleFluidMass(real mass, Chars unit) -> void;
 
    auto scaleSolidMass(real mass, Chars unit) -> void;
 
    // --------------------------------------------------------------------------------------------
    // METHODS TO SCALE THE VOLUME OF THE SYSTEM OR PART OF IT
    // --------------------------------------------------------------------------------------------
 
    auto scaleVolume(real volume, Chars unit) -> void;
 
    auto scalePhaseVolume(StringOrIndex const& phase, real volume, Chars unit) -> void;
 
    auto scaleFluidVolume(real volume, Chars unit) -> void;
 
    auto scaleSolidVolume(real volume, Chars unit) -> void;
 
    // --------------------------------------------------------------------------------------------
    // METHODS FOR UPDATING CHEMICAL STATE AND ITS PROPERTIES
    // --------------------------------------------------------------------------------------------
 
    auto update(real const& T, real const& P, ArrayXrConstRef const& n) -> void;
 
    auto updateIdeal(real const& T, real const& P, ArrayXrConstRef const& n) -> void;
 
    // --------------------------------------------------------------------------------------------
    // MISCELLANEOUS METHODS
    // --------------------------------------------------------------------------------------------
 
    auto system() const -> ChemicalSystem const&;
 
    auto props() const -> ChemicalProps const&;
 
    auto props() -> ChemicalProps&;
 
    auto equilibrium() const -> Equilibrium const&;
 
    auto equilibrium() -> Equilibrium&;
 
    auto output(std::ostream& out) const -> void;
 
    auto output(String const& filename) const -> void;
 
private:
    struct Impl;
 
    Ptr<Impl> pimpl;
};
 
//=================================================================================================
//
// ChemicalState::Equilibrium
//
//=================================================================================================
 
class ChemicalState::Equilibrium
{
public:
    Equilibrium(ChemicalSystem const& system);
 
    Equilibrium(Equilibrium const& other);
 
    virtual ~Equilibrium();
 
    auto operator=(Equilibrium other) -> Equilibrium&;
 
    auto reset() -> void;
 
    auto setNamesInputVariables(Strings const& wnames) -> void;
 
    auto setNamesControlVariablesP(Strings const& pnames) -> void;
 
    auto setNamesControlVariablesQ(Strings const& qnames) -> void;
 
    auto setInputVariables(ArrayXdConstRef const& w) -> void;
 
    auto setControlVariablesP(ArrayXdConstRef const& p) -> void;
 
    auto setControlVariablesQ(ArrayXdConstRef const& q) -> void;
 
    auto setInitialComponentAmounts(ArrayXdConstRef const& c0) -> void;
 
    auto setOptimaState(Optima::State const& state) -> void;
 
    auto empty() const -> bool;
 
    auto numPrimarySpecies() const -> Index;
 
    auto numSecondarySpecies() const -> Index;
 
    auto indicesPrimarySpecies() const -> ArrayXlConstRef;
 
    auto indicesSecondarySpecies() const -> ArrayXlConstRef;
 
    auto elementChemicalPotentials() const -> ArrayXdConstRef;
 
    auto speciesStabilities() const -> ArrayXdConstRef;
 
    auto explicitTitrantAmount(String const& name) const -> real;
 
    auto implicitTitrantAmount(String const& name) const -> real;
 
    auto titrantAmount(String const& name) const -> real;
 
    auto namesInputVariables() const -> Strings const&;
 
    auto namesControlVariablesP() const -> Strings const&;
 
    auto namesControlVariablesQ() const -> Strings const&;
 
    auto inputVariables() const -> ArrayXdConstRef;
 
    auto controlVariablesP() const -> ArrayXdConstRef;
 
    auto controlVariablesQ() const -> ArrayXdConstRef;
 
    auto initialComponentAmounts() const -> ArrayXdConstRef;
 
    auto w() const -> ArrayXdConstRef;
 
    auto p() const -> ArrayXdConstRef;
 
    auto q() const -> ArrayXdConstRef;
 
    auto c() const -> ArrayXdConstRef;
 
    auto optimaState() const -> Optima::State const&;
 
private:
    struct Impl;
 
    Ptr<Impl> pimpl;
};
 
auto operator<<(std::ostream& out, ChemicalState const& state) -> std::ostream&;
 
} // namespace Reaktoro
 
//=========================================================================
// CODE BELOW NEEDED FOR MEMOIZATION TECHNIQUE INVOLVING CHEMICALSTATE
//=========================================================================
 
namespace Reaktoro {
 
template<typename T>
struct MemoizationTraits;
 
template<>
struct MemoizationTraits<ChemicalState>
{
    using Type = ChemicalState;
 
    using CacheType = Tuple<real, real, ArrayXr>;
 
    auto equal(const Tuple<real, real, ArrayXr>& a, ChemicalState const& b)
    {
        auto const& [T, P, n] = a;
        return T == b.temperature() && P == b.pressure() && (n == b.speciesAmounts()).all();
    }
 
    auto assign(Tuple<real, real, ArrayXr>& a, ChemicalState const& b)
    {
        auto& [T, P, n] = a;
        T = b.temperature();
        P = b.pressure();
        n = b.speciesAmounts();
    }
};
 
} // namespace Reaktoro