#![cfg(feature = "knitro_rs")]
use anyhow::Result;
use approx::assert_relative_eq;
use feos::core::parameter::{IdentifierOption, Parameter, PureRecord};
use feos::core::si::{KELVIN, MOL};
use feos::core::{EosResult, Residual, State};
use feos::pcsaft::{PcSaft, PcSaftBinaryRecord, PcSaftParameters, PcSaftRecord};
use feos_campd::process::ProcessModel;
use feos_campd::{
CoMTCAMD, OptimizationProblem, OuterApproximationAlgorithm, PcSaftPropertyModel,
ProcessVariables, Variable,
};
use indexmap::IndexMap;
use ndarray::arr1;
use std::sync::Arc;
#[test]
pub fn test_binary() -> Result<()> {
let camd = CoMTCAMD::from_molecules(vec![
"pentane".into(),
"benzene".into(),
"diethyl ether".into(),
]);
let pcsaft = PcSaftPropertyModel::from_json_molecules(
"tests/gross2001.json",
"tests/poling2000.json",
Some("tests/mixture_test_comps_binary.json"),
IdentifierOption::Name,
true,
)?;
let y0 = camd.get_initial_values("molecule", &IndexMap::from([("pentane", 1)]));
let y1 = camd.get_initial_values("molecule", &IndexMap::from([("diethyl ether", 1)]));
let y = [y0, y1];
let camd_binary = [camd.clone(), camd];
let mut problem = OptimizationProblem::new(camd_binary, pcsaft, CriticalPointModel);
let result = problem.solve_outer_approximation(
y,
OuterApproximationAlgorithm::DuranGrossmann(true),
Some("tests/options_target.opt"),
Some("tests/options_MILP.opt"),
)?;
println!("{} {}", result.target, result.x[0]);
let pentane = PureRecord::new(
Default::default(),
0.0,
PcSaftRecord::new(
2.6896, 3.7729, 231.2, None, None, None, None, None, None, None, None, None, None,
),
);
let benzene = PureRecord::new(
Default::default(),
0.0,
PcSaftRecord::new(
2.4653, 3.6478, 287.35, None, None, None, None, None, None, None, None, None, None,
),
);
let params = PcSaftParameters::new_binary(
vec![benzene, pentane],
Some(PcSaftBinaryRecord::new(
Some(-0.18371738189391876),
None,
None,
)),
)?;
let eos = Arc::new(PcSaft::new(Arc::new(params)));
let (t, _, _) = CriticalPointModel.solve(&eos, &result.x)?;
println!("{t}");
assert_relative_eq!(t, result.target, max_relative = 1e-5);
Ok(())
}
struct CriticalPointModel;
impl<E: Residual> ProcessModel<E> for CriticalPointModel {
fn variables(&self) -> ProcessVariables {
vec![("x", Variable::continuous(0.0, 1.0, 0.1))].into()
}
fn equality_constraints(&self) -> usize {
0
}
fn inequality_constraints(&self) -> usize {
0
}
fn solve(&self, eos: &Arc<E>, x: &[f64]) -> EosResult<(f64, Vec<f64>, Vec<f64>)> {
println!("{x:?}");
let molefracs = arr1(&[x[0], 1.0 - x[0]]);
let cp = State::critical_point(eos, Some(&(molefracs * MOL)), None, Default::default())?;
Ok((-cp.temperature.convert_into(KELVIN), vec![], vec![]))
}
}