Skip to content
Snippets Groups Projects
multifluid_mutant.hpp 5.74 KiB
Newer Older
  • Learn to ignore specific revisions
  • #pragma once
    
    
    namespace teqp {
    
        /**
        This class holds a lightweight reference to the core parts of the model
    
        The reducing and departure functions are moved into this class, while the donor class is used for the corresponding states portion
        */
        template<typename DepartureFunction, typename BaseClass>
        class MultiFluidAdapter {
    
        private:
            std::string meta = "";
    
        public:
            const BaseClass& base;
            const ReducingFunctions redfunc;
    
            const DepartureFunction dep;
    
    
            template<class VecType>
            auto R(const VecType& molefrac) const { return base.R(molefrac); }
    
    
            MultiFluidAdapter(const BaseClass& base, ReducingFunctions&& redfunc, DepartureFunction&& depfunc) : base(base), redfunc(redfunc), dep(depfunc) {};
    
    
            /// Store some sort of metadata in string form (perhaps a JSON representation of the model?)
            void set_meta(const std::string& m) { meta = m; }
            /// Get the metadata stored in string form
            auto get_meta() const { return meta; }
    
            template<typename TType, typename RhoType, typename MoleFracType>
            auto alphar(const TType& T,
                const RhoType& rho,
                const MoleFracType& molefrac) const
            {
                auto Tred = forceeval(redfunc.get_Tr(molefrac));
                auto rhored = forceeval(redfunc.get_rhor(molefrac));
                auto delta = forceeval(rho / rhored);
                auto tau = forceeval(Tred / T);
    
                auto val = base.corr.alphar(tau, delta, molefrac) + dep.alphar(tau, delta, molefrac);
    
                return forceeval(val);
            }
        };
    
        template<class Model>
    
        auto build_multifluid_mutant(const Model& model, const nlohmann::json& jj) {
    
    
            auto N = model.redfunc.Tc.size();
    
            // Allocate the matrices of default models and F factors
            Eigen::MatrixXd F(N, N); F.setZero();
            std::vector<std::vector<DepartureTerms>> funcs(N);
            for (auto i = 0; i < N; ++i) { funcs[i].resize(N); }
    
            // Build the F and departure function matrix
            for (auto i = 0; i < N; ++i) {
                for (auto j = i; j < N; ++j) {
                    if (i == j) {
                        funcs[i][i].add_term(NullEOSTerm());
                    }
                    else {
                        // Extract the given entry
                        auto entry = jj[std::to_string(i)][std::to_string(j)];
    
                        // Set the entry in the matrix of F and departure functions
                        auto dep = entry.at("departure");
                        auto BIP = entry.at("BIP");
                        F(i, j) = BIP.at("Fij");
                        F(j, i) = F(i, j);
                        funcs[i][j] = build_departure_function(dep);
                        funcs[j][i] = build_departure_function(dep);
                    }
                }
            }
    
            // Determine what sort of reducing function is to be used
            auto get_reducing = [&](const auto& deptype) {
                auto red = model.redfunc;
                auto Tc = red.Tc, vc = red.vc;
                if (deptype == "invariant") {
                    using mat = std::decay_t<decltype(MultiFluidInvariantReducingFunction::phiT)>;
                    mat phiT = mat::Zero(N, N), lambdaT = mat::Zero(N, N), phiV = mat::Zero(N, N), lambdaV = mat::Zero(N, N);
    
                    for (auto i = 0; i < N; ++i) {
    
                        for (auto j = i+1; j < N; ++j) {
    
                            // Extract the given entry
    
                            auto entry = jj.at(std::to_string(i)).at(std::to_string(j));
    
    
                            auto BIP = entry.at("BIP");
                            // Set the reducing function parameters in the copy
                            phiT(i, j) = BIP.at("phiT");
                            phiT(j, i) = phiT(i, j);
                            lambdaT(i, j) = BIP.at("lambdaT");
                            lambdaT(j, i) = -lambdaT(i, j);
    
                            phiV(i, j) = BIP.at("phiV");
                            phiV(j, i) = phiV(i, j);
                            lambdaV(i, j) = BIP.at("lambdaV");
                            lambdaV(j, i) = -lambdaV(i, j);
                        }
                    }
                    return ReducingFunctions(MultiFluidInvariantReducingFunction(phiT, lambdaT, phiV, lambdaV, Tc, vc));
                }
                else {
                    using mat = std::decay_t<decltype(MultiFluidReducingFunction::betaT)>;
                    mat betaT = mat::Zero(N, N), gammaT = mat::Zero(N, N), betaV = mat::Zero(N, N), gammaV = mat::Zero(N, N);
    
                    for (auto i = 0; i < N; ++i) {
    
                        for (auto j = i+1; j < N; ++j) {
    
                            // Extract the given entry
    
                            auto entry = jj.at(std::to_string(i)).at(std::to_string(j));
    
                            auto BIP = entry.at("BIP");
                            // Set the reducing function parameters in the copy
                            betaT(i, j) = BIP.at("betaT");
                            betaT(j, i) = 1 / betaT(i, j);
                            betaV(i, j) = BIP.at("betaV");
                            betaV(j, i) = 1 / betaV(i, j);
                            gammaT(i, j) = BIP.at("gammaT"); gammaT(j, i) = gammaT(i, j);
                            gammaV(i, j) = BIP.at("gammaV"); gammaV(j, i) = gammaV(i, j);
                        }
                    }
                    return ReducingFunctions(MultiFluidReducingFunction(betaT, gammaT, betaV, gammaV, Tc, vc));
                }
            };
    
            std::string deptype = (jj.at("0").at("1").at("BIP").contains("type")) ? jj.at("0").at("1").at("BIP")["type"] : "";
    
            ReducingFunctions newred = get_reducing(deptype);
    
    
            auto newdep = DepartureContribution(std::move(F), std::move(funcs));
            auto mfa = MultiFluidAdapter(model, std::move(newred), std::move(newdep));
            /// Store the model spec in the adapted multifluid class
            mfa.set_meta(jj.dump());
            return mfa;
        }