C++ Integration

This example demonstrates how to use fluids from cpp.

Source Code

  1#include <pybind11/pybind11.h>
  2#include <pybind11/embed.h> // Added this for interpreter management
  3#include <pybind11/stl.h>
  4#include <chrono>
  5#include <iostream>
  6#include <iomanip>
  7#include <vector>
  8
  9namespace py = pybind11;
 10PYBIND11_EMBEDDED_MODULE(dummy, m) {} // Needed for embedding
 11
 12class FluidsTester {
 13private:
 14    py::scoped_interpreter guard; // RAII management of the Python interpreter
 15    py::module_ fluids;
 16    
 17public:
 18    FluidsTester() {
 19        try {
 20            // The interpreter is automatically initialized by scoped_interpreter
 21            fluids = py::module_::import("fluids");
 22            std::cout << "✓ Successfully imported fluids\n";
 23            std::cout << "✓ Fluids version: " << fluids.attr("__version__").cast<std::string>() << "\n";
 24        } catch (const std::exception& e) {
 25            std::cerr << "Error initializing Python: " << e.what() << std::endl;
 26            throw;
 27        }
 28    }
 29
 30    void test_fluids() {
 31        try {
 32            // Test basic Reynolds number calculation
 33            auto Re = fluids.attr("Reynolds")(
 34                py::arg("V") = 2.5,
 35                py::arg("D") = 0.1,
 36                py::arg("rho") = 1000,
 37                py::arg("mu") = 0.001
 38            ).cast<double>();
 39            std::cout << "✓ Reynolds number calculation successful: " << Re << "\n";
 40            assert(Re > 0);
 41
 42            // Test friction factor calculation
 43            auto fd = fluids.attr("friction_factor")(
 44                py::arg("Re") = 1e5,
 45                py::arg("eD") = 0.0001
 46            ).cast<double>();
 47            std::cout << "✓ Friction factor calculation successful: " << fd << "\n";
 48            assert(fd > 0 && fd < 1);
 49
 50            std::cout << "\nAll basic tests completed successfully!\n";
 51        } catch (const std::exception& e) {
 52            std::cerr << "Error in fluids tests: " << e.what() << std::endl;
 53            throw;
 54        }
 55    }
 56
 57    void test_atmosphere() {
 58        try {
 59            // Test ATMOSPHERE_1976 class
 60            auto atm = fluids.attr("ATMOSPHERE_1976")(py::arg("Z") = 5000);
 61            
 62            std::cout << "\nTesting atmosphere at 5000m elevation:\n";
 63            std::cout << "✓ Temperature: " << std::fixed << std::setprecision(4) 
 64                      << atm.attr("T").cast<double>() << "\n";
 65            std::cout << "✓ Pressure: " << atm.attr("P").cast<double>() << "\n";
 66            std::cout << "✓ Density: " << std::setprecision(6) 
 67                      << atm.attr("rho").cast<double>() << "\n";
 68            
 69            // Test derived properties
 70            std::cout << "✓ Gravity: " << atm.attr("g").cast<double>() << "\n";
 71            std::cout << "✓ Viscosity: " << std::scientific 
 72                      << atm.attr("mu").cast<double>() << "\n";
 73            std::cout << "✓ Thermal conductivity: " << std::fixed 
 74                      << atm.attr("k").cast<double>() << "\n";
 75            std::cout << "✓ Sonic velocity: " << std::setprecision(4) 
 76                      << atm.attr("v_sonic").cast<double>() << "\n";
 77
 78            // Test static methods
 79            auto atm_class = fluids.attr("ATMOSPHERE_1976");
 80            auto g_high = atm_class.attr("gravity")(py::arg("Z") = 1E5).cast<double>();
 81            std::cout << "✓ High altitude gravity: " << std::setprecision(6) 
 82                      << g_high << "\n";
 83
 84        } catch (const std::exception& e) {
 85            std::cerr << "Error in atmosphere tests: " << e.what() << std::endl;
 86            throw;
 87        }
 88    }
 89
 90    void test_tank() {
 91        try {
 92            // Test basic tank creation
 93            auto T1 = fluids.attr("TANK")(
 94                py::arg("V") = 10,
 95                py::arg("L_over_D") = 0.7,
 96                py::arg("sideB") = "conical",
 97                py::arg("horizontal") = false
 98            );
 99            
100            std::cout << "\nTesting tank calculations:\n";
101            std::cout << "✓ Tank length: " << std::fixed << std::setprecision(6) 
102                      << T1.attr("L").cast<double>() << "\n";
103            std::cout << "✓ Tank diameter: " << T1.attr("D").cast<double>() << "\n";
104
105            // Test torispherical tank
106            auto DIN = fluids.attr("TANK")(
107                py::arg("L") = 3,
108                py::arg("D") = 5,
109                py::arg("horizontal") = false,
110                py::arg("sideA") = "torispherical",
111                py::arg("sideB") = "torispherical",
112                py::arg("sideA_f") = 1,
113                py::arg("sideA_k") = 0.1,
114                py::arg("sideB_f") = 1,
115                py::arg("sideB_k") = 0.1
116            );
117
118            std::cout << "✓ Tank max height: " << DIN.attr("h_max").cast<double>() << "\n";
119            std::cout << "✓ Height at V=40: " << DIN.attr("h_from_V")(40).cast<double>() << "\n";
120            std::cout << "✓ Volume at h=4.1: " << std::setprecision(5) 
121                      << DIN.attr("V_from_h")(4.1).cast<double>() << "\n";
122            
123        } catch (const std::exception& e) {
124            std::cerr << "Error in tank tests: " << e.what() << std::endl;
125            throw;
126        }
127    }
128void test_reynolds() {
129    try {
130        std::cout << "\nTesting Reynolds number calculations:\n";
131        
132        // Test with density and viscosity
133        auto Re1 = fluids.attr("Reynolds")(
134            py::arg("V") = 2.5,
135            py::arg("D") = 0.25,
136            py::arg("rho") = 1.1613,
137            py::arg("mu") = 1.9E-5
138        ).cast<double>();
139        std::cout << "✓ Re (with rho, mu): " << std::fixed << std::setprecision(4) << Re1 << "\n";
140        assert(std::abs(Re1 - 38200.6579) < 0.1);
141        
142        // Test with kinematic viscosity
143        auto Re2 = fluids.attr("Reynolds")(
144            py::arg("V") = 2.5,
145            py::arg("D") = 0.25,
146            py::arg("nu") = 1.636e-05
147        ).cast<double>();
148        std::cout << "✓ Re (with nu): " << Re2 << "\n";
149        assert(std::abs(Re2 - 38202.934) < 0.1);
150        
151    } catch (const std::exception& e) {
152        std::cerr << "Error in Reynolds tests: " << e.what() << std::endl;
153        throw;
154    }
155}
156
157void test_psd() {
158    try {
159        std::cout << "\nTesting particle size distribution functionality:\n";
160        
161        // Create vectors for discrete PSD
162        std::vector<double> ds = {240, 360, 450, 562.5, 703, 878, 1097, 1371, 1713, 
163                                 2141, 2676, 3345, 4181, 5226, 6532};
164        std::vector<double> numbers = {65, 119, 232, 410, 629, 849, 990, 981, 825, 
165                                     579, 297, 111, 21, 1};
166        
167        // Create a discrete PSD
168        auto psd = fluids.attr("particle_size_distribution").attr("ParticleSizeDistribution")(
169            py::arg("ds") = ds,
170            py::arg("fractions") = numbers,
171            py::arg("order") = 0
172        );
173        std::cout << "✓ Created discrete PSD\n";
174        
175        // Test mean sizes
176        auto d21 = psd.attr("mean_size")(2, 1).cast<double>();
177        std::cout << "✓ Size-weighted mean diameter: " << std::fixed 
178                  << std::setprecision(4) << d21 << "\n";
179        assert(std::abs(d21 - 1857.788) < 0.1);
180        
181        auto d10 = psd.attr("mean_size")(1, 0).cast<double>();
182        std::cout << "✓ Arithmetic mean diameter: " << d10 << "\n";
183        assert(std::abs(d10 - 1459.372) < 0.1);
184        
185        // Test percentile calculations
186        auto d10_percentile = psd.attr("dn")(0.1).cast<double>();
187        auto d90_percentile = psd.attr("dn")(0.9).cast<double>();
188        std::cout << "✓ D10: " << d10_percentile << "\n";
189        std::cout << "✓ D90: " << d90_percentile << "\n";
190        
191        // Test probability functions
192        auto pdf_val = psd.attr("pdf")(1000).cast<double>();
193        auto cdf_val = psd.attr("cdf")(5000).cast<double>();
194        std::cout << "✓ PDF at 1000: " << std::scientific << pdf_val << "\n";
195        std::cout << "✓ CDF at 5000: " << std::fixed << std::setprecision(6) << cdf_val << "\n";
196        
197        // Test lognormal distribution
198        auto psd_log = fluids.attr("particle_size_distribution").attr("PSDLognormal")(
199            py::arg("s") = 0.5,
200            py::arg("d_characteristic") = 5E-6
201        );
202        std::cout << "✓ Created lognormal PSD\n";
203        
204        auto vssa = psd_log.attr("vssa").cast<double>();
205        std::cout << "✓ Volume specific surface area: " << std::setprecision(2) 
206                  << vssa << "\n";
207        
208        auto span = psd_log.attr("dn")(0.9).cast<double>() - 
209                   psd_log.attr("dn")(0.1).cast<double>();
210        std::cout << "✓ Span: " << std::scientific << span << "\n";
211        
212        auto ratio_7525 = (psd_log.attr("dn")(0.75).cast<double>() / 
213                          psd_log.attr("dn")(0.25).cast<double>());
214        std::cout << "✓ D75/D25 ratio: " << std::fixed << std::setprecision(6) 
215                  << ratio_7525 << "\n";
216        
217    } catch (const std::exception& e) {
218        std::cerr << "Error in PSD tests: " << e.what() << std::endl;
219        throw;
220    }
221}
222    void benchmark_fluids() {
223        std::cout << "\nRunning benchmarks:\n";
224        
225        // Benchmark friction factor calculation
226        std::cout << "\nBenchmarking friction_factor:\n";
227        auto start = std::chrono::high_resolution_clock::now();
228        
229        for(int i = 0; i < 1000000; i++) {
230            fluids.attr("friction_factor")(
231                py::arg("Re") = 1e5,
232                py::arg("eD") = 0.0001
233            );
234        }
235        
236        auto end = std::chrono::high_resolution_clock::now();
237        auto duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start);
238        
239        std::cout << "Time for 1e6 friction_factor calls: " 
240                  << duration.count() << " microseconds\n";
241        std::cout << "Average time per call: " 
242                  << duration.count() / 1000000.0 << " microseconds\n";
243
244        // Benchmark tank creation
245        std::cout << "\nBenchmarking TANK creation:\n";
246        start = std::chrono::high_resolution_clock::now();
247        
248        for(int i = 0; i < 1000; i++) {
249            fluids.attr("TANK")(
250                py::arg("L") = 3,
251                py::arg("D") = 5,
252                py::arg("horizontal") = false,
253                py::arg("sideA") = "torispherical",
254                py::arg("sideB") = "torispherical",
255                py::arg("sideA_f") = 1,
256                py::arg("sideA_k") = 0.1,
257                py::arg("sideB_f") = 1,
258                py::arg("sideB_k") = 0.1
259            );
260        }
261        
262        end = std::chrono::high_resolution_clock::now();
263        duration = std::chrono::duration_cast<std::chrono::microseconds>(end - start);
264        
265        std::cout << "Average time per tank creation: " 
266                  << duration.count() / 1000.0 << " microseconds\n";
267    }
268};
269
270int main() {
271    try {
272        std::cout << "Running fluids tests from C++...\n";
273        FluidsTester tester;
274        tester.test_fluids();
275        tester.test_atmosphere();
276        tester.test_tank();
277        tester.test_reynolds();
278        tester.test_psd();
279        tester.benchmark_fluids();
280        std::cout << "\nAll tests completed!\n";
281        return 0;
282    } catch (const std::exception& e) {
283        std::cerr << "Fatal error: " << e.what() << std::endl;
284        return 1;
285    }
286}

Requirements

  • Python with fluids installed

  • pybind11

Usage Notes

  • The example demonstrates basic integration with fluids

  • Essentially zero overhead - 1.7 microseconds for friction factor, 8 microseconds for tank creation; only C using the has been observed to have a faster interface, able to save < 0.5 microseconds/call