From 6d663d7a7bb5c25b2bebf20cef233d8cba004504 Mon Sep 17 00:00:00 2001 From: jgostick Date: Sat, 4 Nov 2023 01:33:13 +0900 Subject: [PATCH] comparing conduits to analytical solution in integration test --- tests/integration/Flow_in_straight_conduit.py | 20 ++++++++++++------- 1 file changed, 13 insertions(+), 7 deletions(-) diff --git a/tests/integration/Flow_in_straight_conduit.py b/tests/integration/Flow_in_straight_conduit.py index 8b3473b65a..8042049ab1 100644 --- a/tests/integration/Flow_in_straight_conduit.py +++ b/tests/integration/Flow_in_straight_conduit.py @@ -5,8 +5,13 @@ import openpnm as op import numpy as np + mu = 1.0 + Pin = 1.0 + Pout = 0.0 diameters = [0.9, 1.0, 1.1, 1.5] - for D in diameters: + Q = np.pi*(np.array(diameters)/2)**4/(8*mu*4) * (Pin - Pout) + + for i, D in enumerate(diameters): pn5 = op.network.Cubic(shape=[5, 1, 1], spacing=1) pn5['pore.diameter'] = D pn5['throat.diameter'] = D @@ -24,26 +29,27 @@ pn3.regenerate_models() ph1 = op.phase.Phase(network=pn3) - ph1['pore.viscosity'] = 1.0 + ph1['pore.viscosity'] = mu ph1.add_model( propname='throat.hydraulic_conductance', model=op.models.physics.hydraulic_conductance.generic_hydraulic) ph2 = op.phase.Phase(network=pn5) - ph2['pore.viscosity'] = 1.0 + ph2['pore.viscosity'] = mu ph2.add_model( propname='throat.hydraulic_conductance', model=op.models.physics.hydraulic_conductance.generic_hydraulic) sf1 = op.algorithms.StokesFlow(network=pn3, phase=ph1) - sf1.set_value_BC(pores=pn3.pores('xmin'), values=1.0) - sf1.set_value_BC(pores=pn3.pores('xmax'), values=0.0) + sf1.set_value_BC(pores=pn3.pores('xmin'), values=Pin) + sf1.set_value_BC(pores=pn3.pores('xmax'), values=Pout) sf1.run() r1 = sf1.rate(pores=pn3.pores('xmin')) sf2 = op.algorithms.StokesFlow(network=pn5, phase=ph1) - sf2.set_value_BC(pores=pn5.pores('xmin'), values=1.0) - sf2.set_value_BC(pores=pn5.pores('xmax'), values=0.0) + sf2.set_value_BC(pores=pn5.pores('xmin'), values=Pin) + sf2.set_value_BC(pores=pn5.pores('xmax'), values=Pout) sf2.run() r2 = sf2.rate(pores=pn5.pores('xmin')) assert np.allclose(r1, r2, atol=0, rtol=1e-10) + assert np.allclose(r1, Q[i], atol=0, rtol=1e-10)