From c977d5439960263bc7c535a611cc8afc4d4a0392 Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Fri, 13 Mar 2026 09:58:43 +0100 Subject: [PATCH 1/8] Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py --- avaframe/com1DFA/DFAfunctionsCython.pyx | 9 +++----- avaframe/com1DFA/checkCfg.py | 3 +-- avaframe/com1DFA/com1DFA.py | 3 +-- avaframe/com1DFA/com1DFACfg.ini | 8 +++----- avaframe/tests/test_DFAfunctionsCython.py | 8 ++------ docs/theoryCom1DFA.rst | 25 +++++++++++++++-------- 6 files changed, 27 insertions(+), 29 deletions(-) diff --git a/avaframe/com1DFA/DFAfunctionsCython.pyx b/avaframe/com1DFA/DFAfunctionsCython.pyx index 72b859b23..3ef93f80b 100644 --- a/avaframe/com1DFA/DFAfunctionsCython.pyx +++ b/avaframe/com1DFA/DFAfunctionsCython.pyx @@ -95,10 +95,9 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double alpha2TauyObrienAndJulien = cfg.getfloat('alpha2TauyObrienAndJulien') cdef double beta2TauyObrienAndJulien = cfg.getfloat('beta2TauyObrienAndJulien') cdef double alphaObrienAndJulien = cfg.getfloat('alphaObrienAndJulien') - cdef double alpha1EtaHerschelAndBulkley = cfg.getfloat('alpha1EtaHerschelAndBulkley') - cdef double beta1EtaHerschelAndBulkley = cfg.getfloat('beta1EtaHerschelAndBulkley') cdef double alpha2TauyHerschelAndBulkley = cfg.getfloat('alpha2TauyHerschelAndBulkley') cdef double beta2TauyHerschelAndBulkley = cfg.getfloat('beta2TauyHerschelAndBulkley') + cdef double kHerschelAndBulkley = cfg.getfloat('kHerschelAndBulkley') cdef double nHerschelAndBulkley = cfg.getfloat('nHerschelAndBulkley') cdef double alpha1EtaBingham = cfg.getfloat('alpha1EtaBingham') cdef double beta1EtaBingham = cfg.getfloat('beta1EtaBingham') @@ -170,7 +169,7 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double gravAccNorm, accNormCurv, effAccNorm, gravAccTangX, gravAccTangY, gravAccTangZ, forceBotTang, sigmaB, tau cdef double muVoellmyRaster, xsiVoellmyRaster cdef double shearRate, etaObrienAndJulien, tauyObrienAndJulien, lmObrienAndJulien - cdef double lambdaBagnold, cObrienAndJulien, etaHerschelAndBulkley, tauyHerschelAndBulkley, etaBingham, tauyBingham + cdef double lambdaBagnold, cObrienAndJulien, tauyHerschelAndBulkley, etaBingham, tauyBingham # variables for interpolation cdef int Lx0, Ly0, LxEnd0, LyEnd0, iCell, iCellEnd cdef double w[4] @@ -340,12 +339,10 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType tau = tauyObrienAndJulien + etaObrienAndJulien * shearRate + cObrienAndJulien * (shearRate * shearRate) elif frictType == 11: ## Herschel and Bulkley - # viscosity - etaHerschelAndBulkley = alpha1EtaHerschelAndBulkley * math.exp(beta1EtaHerschelAndBulkley * cvSediment) # yield shear stress tauyHerschelAndBulkley = alpha2TauyHerschelAndBulkley * math.exp(beta2TauyHerschelAndBulkley * cvSediment) # shear stress - tau = tauyHerschelAndBulkley + etaHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) + tau = tauyHerschelAndBulkley + kHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) elif frictType == 12: ## Bingham # viscosity diff --git a/avaframe/com1DFA/checkCfg.py b/avaframe/com1DFA/checkCfg.py index ef8289a43..83c76f780 100644 --- a/avaframe/com1DFA/checkCfg.py +++ b/avaframe/com1DFA/checkCfg.py @@ -135,10 +135,9 @@ def checkCfgFrictionModel(cfg, inputSimFiles, relVolume=""): "alpha2TauyObrienAndJulien", "beta2TauyObrienAndJulien", "alphaObrienAndJulien", - "alpha1EtaHerschelAndBulkley", - "beta1EtaHerschelAndBulkley", "alpha2TauyHerschelAndBulkley", "beta2TauyHerschelAndBulkley", + "kHerschelAndBulkley", "nHerschelAndBulkley", "alpha1EtaBingham", "beta1EtaBingham", diff --git a/avaframe/com1DFA/com1DFA.py b/avaframe/com1DFA/com1DFA.py index ef2993f9a..20cad75a2 100644 --- a/avaframe/com1DFA/com1DFA.py +++ b/avaframe/com1DFA/com1DFA.py @@ -932,8 +932,7 @@ def createReportDict(avaDir, logName, relName, inputSimLines, cfg, reportAreaInf "type": "columns", "model": "Herschel and Bulkley", "n": cfgGen["nHerschelAndBulkley"], - "alpha1": cfgGen["alpha1EtaHerschelAndBulkley"], - "beta1": cfgGen["beta1EtaHerschelAndBulkley"], + "K": cfgGen["kHerschelAndBulkley"], "alpha2": cfgGen["alpha2TauyHerschelAndBulkley"], "beta2": cfgGen["beta2TauyHerschelAndBulkley"], } diff --git a/avaframe/com1DFA/com1DFACfg.ini b/avaframe/com1DFA/com1DFACfg.ini index 1a47466e1..dae998349 100644 --- a/avaframe/com1DFA/com1DFACfg.ini +++ b/avaframe/com1DFA/com1DFACfg.ini @@ -401,16 +401,14 @@ beta2TauyObrienAndJulien = 13.11 # parameter, default value Takahashi (1978) alphaObrienAndJulien = 0.01 #++++++++++++Herschel and Bulkley -# viscosity parameter -## TODO adjusting, default value O´Brien and Julien (1985) -alpha1EtaHerschelAndBulkley = 0.650 - ## TODO adjusting, default value O´Brien and Julien (1985) -beta1EtaHerschelAndBulkley = 16.81 # yield shear stress parameter ## TODO adjusting, default value O´Brien and Julien (1985) alpha2TauyHerschelAndBulkley = 0.00886 ## TODO adjusting, default value O´Brien and Julien (1985) beta2TauyHerschelAndBulkley = 13.11 +# consistency factor K +## TODO adjusting +kHerschelAndBulkley = 50 # flow exponent ## TODO adjusting nHerschelAndBulkley = 1.5 diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 816ce665b..488ee2c05 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -985,19 +985,15 @@ def test_computeForceC(): ## Herschel and Bulkley # read input parameters - test_alpha1EtaHerschelAndBulkley = float(cfg["GENERAL"]["alpha1EtaHerschelAndBulkley"]) - test_beta1EtaHerschelAndBulkley = float(cfg["GENERAL"]["beta1EtaHerschelAndBulkley"]) test_alpha2TauyHerschelAndBulkley = float(cfg["GENERAL"]["alpha2TauyHerschelAndBulkley"]) test_beta2TauyHerschelAndBulkley = float(cfg["GENERAL"]["beta2TauyHerschelAndBulkley"]) + test_kHerschelAndBulkley = float(cfg["GENERAL"]["kHerschelAndBulkley"]) test_nHerschelAndBulkley = float(cfg["GENERAL"]["nHerschelAndBulkley"]) # compute tau - test_etaHerschelandBulkley = test_alpha1EtaHerschelAndBulkley * math.exp( - test_beta1EtaHerschelAndBulkley * test_cvSediment - ) test_tauyHerschelandBulkley = test_alpha2TauyHerschelAndBulkley * math.exp( test_beta2TauyHerschelAndBulkley * test_cvSediment ) - test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_etaHerschelandBulkley * math.pow( + test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_kHerschelAndBulkley * math.pow( test_shearRate[0], test_nHerschelAndBulkley ) test_tauArray[1] = test_tauHerschelAndBulkley diff --git a/docs/theoryCom1DFA.rst b/docs/theoryCom1DFA.rst index da3fce612..ceca2977d 100644 --- a/docs/theoryCom1DFA.rst +++ b/docs/theoryCom1DFA.rst @@ -480,10 +480,11 @@ With &\mathbf{x} \qquad &\text{position vector}\end{aligned} For solid-fluid mixtures like debris flows, where the properties of the fluid phase are dominating the flow process, -the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). +the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` (resp. consistency factor :math:`K` for :ref:`herschelandbulkley`) +and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). .. math:: - \tau^{(b)}_i = f(\eta_m,\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) + \tau^{(b)}_i = f(\eta_m,(K),\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) :label: rheological model With @@ -491,6 +492,7 @@ With .. math:: \begin{aligned} &\eta_m \qquad &\text{dynamic viscosity of the solid-fluid mixture}\\ + &K \qquad &\text{consistency factor of the solid-fluid mixture}\\ &\dot\gamma \qquad &\text{shear rate}\\ &\overline{h} \qquad &\text{average flow thickness}\\ &\rho_m \qquad &\text{density of the solid-fluid mixture}\\ @@ -662,12 +664,13 @@ and non-Newtonian fluids, in which a yield shear stress has to be exceeded or sh The rheological models are incorporated into a single, general form, which can be expressed as follows: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot\gamma^n + C \cdot \dot\gamma^2 + \tau^{(b)} = \tau_y + K \cdot \dot\gamma^n + C \cdot \dot\gamma^2 :label: rheology-general The yield shear stress :math:`\tau_y` :math:`[Pa]` defines a lower limit below which no flow takes place (see :math:`\tau_0` at :ref:`samosatfrict`). -The dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. -:math:`\dot\gamma` is the flow velocity gradient, or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). +:math:`K` is a consistency factor :math:`[Pa \cdot s^{n}]`. For the :ref:`bingham` and the :ref:`obrienandjulien` rheologies the consistency factor corresponds to the dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` +which quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. :math:`\dot\gamma` is the flow velocity gradient, +or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). The flow index :math:`n` describes the rheological behaviour of the mixture as :cite:`KaZhHaHe2025`: - :math:`n = 0` a rate-independent solid-like behaviour, @@ -731,6 +734,8 @@ in following substitution :cite:`Iv1997, DeIv2001, GiFlSaHe2020`: where :math:`\overline{u}` is the depth-averaged flow velocity and :math:`\overline{h}` is the flow thickness. +.. _obrienandjulien: + O´Brien and Julien +++++++++++++++++++ The quadratic rheological model proposed by O´Brien and Julien :cite:`ObJu1985` accounts for various shear stress components, including @@ -761,17 +766,21 @@ of sediment :math:`[kg \cdot m^{-3}]`, :math:`d_s` is the sediment size :math:`[ The quadratic rheological model is particularly suitable for simulating debris flows with high sediment concentrations, in which energy dissipation and resistance due to turbulence and particle collisions are significant. +.. _herschelandbulkley: + Herschel and Bulkley +++++++++++++++++++++++++ The model by Herschel and Bulkley :cite:`HeBu1926` is expressed by an empirical power-law equation: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot{\gamma}^n + \tau^{(b)} = \tau_y + K \cdot \dot{\gamma}^n :label: herschelAndBulkley where the flow index :math:`n` describes the rheological behaviour of the mixture (see above). The factor in front of the shear rate -was originally introduced as a consistency factor :math:`K`. In :py:mod:`com1DFA` :math:`K` equals to the bulk dynamic viscosity :math:`\eta_m`. -This rheology applies to fine-grained soil-water mixtures that exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. +was originally introduced as a consistency factor :math:`K`.This rheology applies to fine-grained soil-water mixtures that +exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. + +.. _bingham: Bingham ++++++++++++++++ From a58bd38b9e20695b05617487c1295ebfd5065807 Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Fri, 13 Mar 2026 09:58:43 +0100 Subject: [PATCH 2/8] Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py --- avaframe/com1DFA/DFAfunctionsCython.pyx | 9 +++----- avaframe/com1DFA/checkCfg.py | 3 +-- avaframe/com1DFA/com1DFA.py | 3 +-- avaframe/com1DFA/com1DFACfg.ini | 8 +++----- avaframe/tests/test_DFAfunctionsCython.py | 8 ++------ docs/theoryCom1DFA.rst | 25 +++++++++++++++-------- 6 files changed, 27 insertions(+), 29 deletions(-) diff --git a/avaframe/com1DFA/DFAfunctionsCython.pyx b/avaframe/com1DFA/DFAfunctionsCython.pyx index 72b859b23..3ef93f80b 100644 --- a/avaframe/com1DFA/DFAfunctionsCython.pyx +++ b/avaframe/com1DFA/DFAfunctionsCython.pyx @@ -95,10 +95,9 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double alpha2TauyObrienAndJulien = cfg.getfloat('alpha2TauyObrienAndJulien') cdef double beta2TauyObrienAndJulien = cfg.getfloat('beta2TauyObrienAndJulien') cdef double alphaObrienAndJulien = cfg.getfloat('alphaObrienAndJulien') - cdef double alpha1EtaHerschelAndBulkley = cfg.getfloat('alpha1EtaHerschelAndBulkley') - cdef double beta1EtaHerschelAndBulkley = cfg.getfloat('beta1EtaHerschelAndBulkley') cdef double alpha2TauyHerschelAndBulkley = cfg.getfloat('alpha2TauyHerschelAndBulkley') cdef double beta2TauyHerschelAndBulkley = cfg.getfloat('beta2TauyHerschelAndBulkley') + cdef double kHerschelAndBulkley = cfg.getfloat('kHerschelAndBulkley') cdef double nHerschelAndBulkley = cfg.getfloat('nHerschelAndBulkley') cdef double alpha1EtaBingham = cfg.getfloat('alpha1EtaBingham') cdef double beta1EtaBingham = cfg.getfloat('beta1EtaBingham') @@ -170,7 +169,7 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double gravAccNorm, accNormCurv, effAccNorm, gravAccTangX, gravAccTangY, gravAccTangZ, forceBotTang, sigmaB, tau cdef double muVoellmyRaster, xsiVoellmyRaster cdef double shearRate, etaObrienAndJulien, tauyObrienAndJulien, lmObrienAndJulien - cdef double lambdaBagnold, cObrienAndJulien, etaHerschelAndBulkley, tauyHerschelAndBulkley, etaBingham, tauyBingham + cdef double lambdaBagnold, cObrienAndJulien, tauyHerschelAndBulkley, etaBingham, tauyBingham # variables for interpolation cdef int Lx0, Ly0, LxEnd0, LyEnd0, iCell, iCellEnd cdef double w[4] @@ -340,12 +339,10 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType tau = tauyObrienAndJulien + etaObrienAndJulien * shearRate + cObrienAndJulien * (shearRate * shearRate) elif frictType == 11: ## Herschel and Bulkley - # viscosity - etaHerschelAndBulkley = alpha1EtaHerschelAndBulkley * math.exp(beta1EtaHerschelAndBulkley * cvSediment) # yield shear stress tauyHerschelAndBulkley = alpha2TauyHerschelAndBulkley * math.exp(beta2TauyHerschelAndBulkley * cvSediment) # shear stress - tau = tauyHerschelAndBulkley + etaHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) + tau = tauyHerschelAndBulkley + kHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) elif frictType == 12: ## Bingham # viscosity diff --git a/avaframe/com1DFA/checkCfg.py b/avaframe/com1DFA/checkCfg.py index ef8289a43..83c76f780 100644 --- a/avaframe/com1DFA/checkCfg.py +++ b/avaframe/com1DFA/checkCfg.py @@ -135,10 +135,9 @@ def checkCfgFrictionModel(cfg, inputSimFiles, relVolume=""): "alpha2TauyObrienAndJulien", "beta2TauyObrienAndJulien", "alphaObrienAndJulien", - "alpha1EtaHerschelAndBulkley", - "beta1EtaHerschelAndBulkley", "alpha2TauyHerschelAndBulkley", "beta2TauyHerschelAndBulkley", + "kHerschelAndBulkley", "nHerschelAndBulkley", "alpha1EtaBingham", "beta1EtaBingham", diff --git a/avaframe/com1DFA/com1DFA.py b/avaframe/com1DFA/com1DFA.py index ef2993f9a..20cad75a2 100644 --- a/avaframe/com1DFA/com1DFA.py +++ b/avaframe/com1DFA/com1DFA.py @@ -932,8 +932,7 @@ def createReportDict(avaDir, logName, relName, inputSimLines, cfg, reportAreaInf "type": "columns", "model": "Herschel and Bulkley", "n": cfgGen["nHerschelAndBulkley"], - "alpha1": cfgGen["alpha1EtaHerschelAndBulkley"], - "beta1": cfgGen["beta1EtaHerschelAndBulkley"], + "K": cfgGen["kHerschelAndBulkley"], "alpha2": cfgGen["alpha2TauyHerschelAndBulkley"], "beta2": cfgGen["beta2TauyHerschelAndBulkley"], } diff --git a/avaframe/com1DFA/com1DFACfg.ini b/avaframe/com1DFA/com1DFACfg.ini index 1a47466e1..dae998349 100644 --- a/avaframe/com1DFA/com1DFACfg.ini +++ b/avaframe/com1DFA/com1DFACfg.ini @@ -401,16 +401,14 @@ beta2TauyObrienAndJulien = 13.11 # parameter, default value Takahashi (1978) alphaObrienAndJulien = 0.01 #++++++++++++Herschel and Bulkley -# viscosity parameter -## TODO adjusting, default value O´Brien and Julien (1985) -alpha1EtaHerschelAndBulkley = 0.650 - ## TODO adjusting, default value O´Brien and Julien (1985) -beta1EtaHerschelAndBulkley = 16.81 # yield shear stress parameter ## TODO adjusting, default value O´Brien and Julien (1985) alpha2TauyHerschelAndBulkley = 0.00886 ## TODO adjusting, default value O´Brien and Julien (1985) beta2TauyHerschelAndBulkley = 13.11 +# consistency factor K +## TODO adjusting +kHerschelAndBulkley = 50 # flow exponent ## TODO adjusting nHerschelAndBulkley = 1.5 diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 816ce665b..488ee2c05 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -985,19 +985,15 @@ def test_computeForceC(): ## Herschel and Bulkley # read input parameters - test_alpha1EtaHerschelAndBulkley = float(cfg["GENERAL"]["alpha1EtaHerschelAndBulkley"]) - test_beta1EtaHerschelAndBulkley = float(cfg["GENERAL"]["beta1EtaHerschelAndBulkley"]) test_alpha2TauyHerschelAndBulkley = float(cfg["GENERAL"]["alpha2TauyHerschelAndBulkley"]) test_beta2TauyHerschelAndBulkley = float(cfg["GENERAL"]["beta2TauyHerschelAndBulkley"]) + test_kHerschelAndBulkley = float(cfg["GENERAL"]["kHerschelAndBulkley"]) test_nHerschelAndBulkley = float(cfg["GENERAL"]["nHerschelAndBulkley"]) # compute tau - test_etaHerschelandBulkley = test_alpha1EtaHerschelAndBulkley * math.exp( - test_beta1EtaHerschelAndBulkley * test_cvSediment - ) test_tauyHerschelandBulkley = test_alpha2TauyHerschelAndBulkley * math.exp( test_beta2TauyHerschelAndBulkley * test_cvSediment ) - test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_etaHerschelandBulkley * math.pow( + test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_kHerschelAndBulkley * math.pow( test_shearRate[0], test_nHerschelAndBulkley ) test_tauArray[1] = test_tauHerschelAndBulkley diff --git a/docs/theoryCom1DFA.rst b/docs/theoryCom1DFA.rst index da3fce612..ceca2977d 100644 --- a/docs/theoryCom1DFA.rst +++ b/docs/theoryCom1DFA.rst @@ -480,10 +480,11 @@ With &\mathbf{x} \qquad &\text{position vector}\end{aligned} For solid-fluid mixtures like debris flows, where the properties of the fluid phase are dominating the flow process, -the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). +the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` (resp. consistency factor :math:`K` for :ref:`herschelandbulkley`) +and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). .. math:: - \tau^{(b)}_i = f(\eta_m,\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) + \tau^{(b)}_i = f(\eta_m,(K),\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) :label: rheological model With @@ -491,6 +492,7 @@ With .. math:: \begin{aligned} &\eta_m \qquad &\text{dynamic viscosity of the solid-fluid mixture}\\ + &K \qquad &\text{consistency factor of the solid-fluid mixture}\\ &\dot\gamma \qquad &\text{shear rate}\\ &\overline{h} \qquad &\text{average flow thickness}\\ &\rho_m \qquad &\text{density of the solid-fluid mixture}\\ @@ -662,12 +664,13 @@ and non-Newtonian fluids, in which a yield shear stress has to be exceeded or sh The rheological models are incorporated into a single, general form, which can be expressed as follows: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot\gamma^n + C \cdot \dot\gamma^2 + \tau^{(b)} = \tau_y + K \cdot \dot\gamma^n + C \cdot \dot\gamma^2 :label: rheology-general The yield shear stress :math:`\tau_y` :math:`[Pa]` defines a lower limit below which no flow takes place (see :math:`\tau_0` at :ref:`samosatfrict`). -The dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. -:math:`\dot\gamma` is the flow velocity gradient, or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). +:math:`K` is a consistency factor :math:`[Pa \cdot s^{n}]`. For the :ref:`bingham` and the :ref:`obrienandjulien` rheologies the consistency factor corresponds to the dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` +which quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. :math:`\dot\gamma` is the flow velocity gradient, +or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). The flow index :math:`n` describes the rheological behaviour of the mixture as :cite:`KaZhHaHe2025`: - :math:`n = 0` a rate-independent solid-like behaviour, @@ -731,6 +734,8 @@ in following substitution :cite:`Iv1997, DeIv2001, GiFlSaHe2020`: where :math:`\overline{u}` is the depth-averaged flow velocity and :math:`\overline{h}` is the flow thickness. +.. _obrienandjulien: + O´Brien and Julien +++++++++++++++++++ The quadratic rheological model proposed by O´Brien and Julien :cite:`ObJu1985` accounts for various shear stress components, including @@ -761,17 +766,21 @@ of sediment :math:`[kg \cdot m^{-3}]`, :math:`d_s` is the sediment size :math:`[ The quadratic rheological model is particularly suitable for simulating debris flows with high sediment concentrations, in which energy dissipation and resistance due to turbulence and particle collisions are significant. +.. _herschelandbulkley: + Herschel and Bulkley +++++++++++++++++++++++++ The model by Herschel and Bulkley :cite:`HeBu1926` is expressed by an empirical power-law equation: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot{\gamma}^n + \tau^{(b)} = \tau_y + K \cdot \dot{\gamma}^n :label: herschelAndBulkley where the flow index :math:`n` describes the rheological behaviour of the mixture (see above). The factor in front of the shear rate -was originally introduced as a consistency factor :math:`K`. In :py:mod:`com1DFA` :math:`K` equals to the bulk dynamic viscosity :math:`\eta_m`. -This rheology applies to fine-grained soil-water mixtures that exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. +was originally introduced as a consistency factor :math:`K`.This rheology applies to fine-grained soil-water mixtures that +exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. + +.. _bingham: Bingham ++++++++++++++++ From b32ed229ef3b2b2872d2c61eff320a35da09cf52 Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Mon, 16 Mar 2026 09:00:16 +0100 Subject: [PATCH 3/8] Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py --- avaframe/tests/test_DFAfunctionsCython.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 488ee2c05..731dec875 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -839,10 +839,9 @@ def test_computeForceC(): "mucoulomb": "0.155", "mucoulombminshear": "0.155", "tau0coulombminshear": "70", - "alpha1EtaHerschelAndBulkley": "0.650", - "beta1EtaHerschelAndBulkley": "16.81", "alpha2TauyHerschelAndBulkley": "0.00886", "beta2TauyHerschelAndBulkley": "13.11", + "kHerschelAndBulkley": "50", "nHerschelAndBulkley": "1.5", "alpha1EtaBingham": "0.650", "beta1EtaBingham": "16.81", From ef010436e07c40085cf6745f3011d7ae0e518bcb Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Mon, 16 Mar 2026 09:00:16 +0100 Subject: [PATCH 4/8] Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py --- avaframe/tests/test_DFAfunctionsCython.py | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 488ee2c05..731dec875 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -839,10 +839,9 @@ def test_computeForceC(): "mucoulomb": "0.155", "mucoulombminshear": "0.155", "tau0coulombminshear": "70", - "alpha1EtaHerschelAndBulkley": "0.650", - "beta1EtaHerschelAndBulkley": "16.81", "alpha2TauyHerschelAndBulkley": "0.00886", "beta2TauyHerschelAndBulkley": "13.11", + "kHerschelAndBulkley": "50", "nHerschelAndBulkley": "1.5", "alpha1EtaBingham": "0.650", "beta1EtaBingham": "16.81", From 3d6e7b4a6e044a6bd3730b8295ec60c4243e49d9 Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Fri, 13 Mar 2026 09:58:43 +0100 Subject: [PATCH 5/8] Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py --- avaframe/com1DFA/DFAfunctionsCython.pyx | 9 +++----- avaframe/com1DFA/checkCfg.py | 3 +-- avaframe/com1DFA/com1DFA.py | 3 +-- avaframe/com1DFA/com1DFACfg.ini | 8 +++----- avaframe/tests/test_DFAfunctionsCython.py | 11 +++------- docs/theoryCom1DFA.rst | 25 +++++++++++++++-------- 6 files changed, 28 insertions(+), 31 deletions(-) diff --git a/avaframe/com1DFA/DFAfunctionsCython.pyx b/avaframe/com1DFA/DFAfunctionsCython.pyx index 72b859b23..3ef93f80b 100644 --- a/avaframe/com1DFA/DFAfunctionsCython.pyx +++ b/avaframe/com1DFA/DFAfunctionsCython.pyx @@ -95,10 +95,9 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double alpha2TauyObrienAndJulien = cfg.getfloat('alpha2TauyObrienAndJulien') cdef double beta2TauyObrienAndJulien = cfg.getfloat('beta2TauyObrienAndJulien') cdef double alphaObrienAndJulien = cfg.getfloat('alphaObrienAndJulien') - cdef double alpha1EtaHerschelAndBulkley = cfg.getfloat('alpha1EtaHerschelAndBulkley') - cdef double beta1EtaHerschelAndBulkley = cfg.getfloat('beta1EtaHerschelAndBulkley') cdef double alpha2TauyHerschelAndBulkley = cfg.getfloat('alpha2TauyHerschelAndBulkley') cdef double beta2TauyHerschelAndBulkley = cfg.getfloat('beta2TauyHerschelAndBulkley') + cdef double kHerschelAndBulkley = cfg.getfloat('kHerschelAndBulkley') cdef double nHerschelAndBulkley = cfg.getfloat('nHerschelAndBulkley') cdef double alpha1EtaBingham = cfg.getfloat('alpha1EtaBingham') cdef double beta1EtaBingham = cfg.getfloat('beta1EtaBingham') @@ -170,7 +169,7 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double gravAccNorm, accNormCurv, effAccNorm, gravAccTangX, gravAccTangY, gravAccTangZ, forceBotTang, sigmaB, tau cdef double muVoellmyRaster, xsiVoellmyRaster cdef double shearRate, etaObrienAndJulien, tauyObrienAndJulien, lmObrienAndJulien - cdef double lambdaBagnold, cObrienAndJulien, etaHerschelAndBulkley, tauyHerschelAndBulkley, etaBingham, tauyBingham + cdef double lambdaBagnold, cObrienAndJulien, tauyHerschelAndBulkley, etaBingham, tauyBingham # variables for interpolation cdef int Lx0, Ly0, LxEnd0, LyEnd0, iCell, iCellEnd cdef double w[4] @@ -340,12 +339,10 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType tau = tauyObrienAndJulien + etaObrienAndJulien * shearRate + cObrienAndJulien * (shearRate * shearRate) elif frictType == 11: ## Herschel and Bulkley - # viscosity - etaHerschelAndBulkley = alpha1EtaHerschelAndBulkley * math.exp(beta1EtaHerschelAndBulkley * cvSediment) # yield shear stress tauyHerschelAndBulkley = alpha2TauyHerschelAndBulkley * math.exp(beta2TauyHerschelAndBulkley * cvSediment) # shear stress - tau = tauyHerschelAndBulkley + etaHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) + tau = tauyHerschelAndBulkley + kHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) elif frictType == 12: ## Bingham # viscosity diff --git a/avaframe/com1DFA/checkCfg.py b/avaframe/com1DFA/checkCfg.py index ef8289a43..83c76f780 100644 --- a/avaframe/com1DFA/checkCfg.py +++ b/avaframe/com1DFA/checkCfg.py @@ -135,10 +135,9 @@ def checkCfgFrictionModel(cfg, inputSimFiles, relVolume=""): "alpha2TauyObrienAndJulien", "beta2TauyObrienAndJulien", "alphaObrienAndJulien", - "alpha1EtaHerschelAndBulkley", - "beta1EtaHerschelAndBulkley", "alpha2TauyHerschelAndBulkley", "beta2TauyHerschelAndBulkley", + "kHerschelAndBulkley", "nHerschelAndBulkley", "alpha1EtaBingham", "beta1EtaBingham", diff --git a/avaframe/com1DFA/com1DFA.py b/avaframe/com1DFA/com1DFA.py index ef2993f9a..20cad75a2 100644 --- a/avaframe/com1DFA/com1DFA.py +++ b/avaframe/com1DFA/com1DFA.py @@ -932,8 +932,7 @@ def createReportDict(avaDir, logName, relName, inputSimLines, cfg, reportAreaInf "type": "columns", "model": "Herschel and Bulkley", "n": cfgGen["nHerschelAndBulkley"], - "alpha1": cfgGen["alpha1EtaHerschelAndBulkley"], - "beta1": cfgGen["beta1EtaHerschelAndBulkley"], + "K": cfgGen["kHerschelAndBulkley"], "alpha2": cfgGen["alpha2TauyHerschelAndBulkley"], "beta2": cfgGen["beta2TauyHerschelAndBulkley"], } diff --git a/avaframe/com1DFA/com1DFACfg.ini b/avaframe/com1DFA/com1DFACfg.ini index 1a47466e1..dae998349 100644 --- a/avaframe/com1DFA/com1DFACfg.ini +++ b/avaframe/com1DFA/com1DFACfg.ini @@ -401,16 +401,14 @@ beta2TauyObrienAndJulien = 13.11 # parameter, default value Takahashi (1978) alphaObrienAndJulien = 0.01 #++++++++++++Herschel and Bulkley -# viscosity parameter -## TODO adjusting, default value O´Brien and Julien (1985) -alpha1EtaHerschelAndBulkley = 0.650 - ## TODO adjusting, default value O´Brien and Julien (1985) -beta1EtaHerschelAndBulkley = 16.81 # yield shear stress parameter ## TODO adjusting, default value O´Brien and Julien (1985) alpha2TauyHerschelAndBulkley = 0.00886 ## TODO adjusting, default value O´Brien and Julien (1985) beta2TauyHerschelAndBulkley = 13.11 +# consistency factor K +## TODO adjusting +kHerschelAndBulkley = 50 # flow exponent ## TODO adjusting nHerschelAndBulkley = 1.5 diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 816ce665b..731dec875 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -839,10 +839,9 @@ def test_computeForceC(): "mucoulomb": "0.155", "mucoulombminshear": "0.155", "tau0coulombminshear": "70", - "alpha1EtaHerschelAndBulkley": "0.650", - "beta1EtaHerschelAndBulkley": "16.81", "alpha2TauyHerschelAndBulkley": "0.00886", "beta2TauyHerschelAndBulkley": "13.11", + "kHerschelAndBulkley": "50", "nHerschelAndBulkley": "1.5", "alpha1EtaBingham": "0.650", "beta1EtaBingham": "16.81", @@ -985,19 +984,15 @@ def test_computeForceC(): ## Herschel and Bulkley # read input parameters - test_alpha1EtaHerschelAndBulkley = float(cfg["GENERAL"]["alpha1EtaHerschelAndBulkley"]) - test_beta1EtaHerschelAndBulkley = float(cfg["GENERAL"]["beta1EtaHerschelAndBulkley"]) test_alpha2TauyHerschelAndBulkley = float(cfg["GENERAL"]["alpha2TauyHerschelAndBulkley"]) test_beta2TauyHerschelAndBulkley = float(cfg["GENERAL"]["beta2TauyHerschelAndBulkley"]) + test_kHerschelAndBulkley = float(cfg["GENERAL"]["kHerschelAndBulkley"]) test_nHerschelAndBulkley = float(cfg["GENERAL"]["nHerschelAndBulkley"]) # compute tau - test_etaHerschelandBulkley = test_alpha1EtaHerschelAndBulkley * math.exp( - test_beta1EtaHerschelAndBulkley * test_cvSediment - ) test_tauyHerschelandBulkley = test_alpha2TauyHerschelAndBulkley * math.exp( test_beta2TauyHerschelAndBulkley * test_cvSediment ) - test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_etaHerschelandBulkley * math.pow( + test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_kHerschelAndBulkley * math.pow( test_shearRate[0], test_nHerschelAndBulkley ) test_tauArray[1] = test_tauHerschelAndBulkley diff --git a/docs/theoryCom1DFA.rst b/docs/theoryCom1DFA.rst index da3fce612..ceca2977d 100644 --- a/docs/theoryCom1DFA.rst +++ b/docs/theoryCom1DFA.rst @@ -480,10 +480,11 @@ With &\mathbf{x} \qquad &\text{position vector}\end{aligned} For solid-fluid mixtures like debris flows, where the properties of the fluid phase are dominating the flow process, -the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). +the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` (resp. consistency factor :math:`K` for :ref:`herschelandbulkley`) +and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). .. math:: - \tau^{(b)}_i = f(\eta_m,\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) + \tau^{(b)}_i = f(\eta_m,(K),\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) :label: rheological model With @@ -491,6 +492,7 @@ With .. math:: \begin{aligned} &\eta_m \qquad &\text{dynamic viscosity of the solid-fluid mixture}\\ + &K \qquad &\text{consistency factor of the solid-fluid mixture}\\ &\dot\gamma \qquad &\text{shear rate}\\ &\overline{h} \qquad &\text{average flow thickness}\\ &\rho_m \qquad &\text{density of the solid-fluid mixture}\\ @@ -662,12 +664,13 @@ and non-Newtonian fluids, in which a yield shear stress has to be exceeded or sh The rheological models are incorporated into a single, general form, which can be expressed as follows: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot\gamma^n + C \cdot \dot\gamma^2 + \tau^{(b)} = \tau_y + K \cdot \dot\gamma^n + C \cdot \dot\gamma^2 :label: rheology-general The yield shear stress :math:`\tau_y` :math:`[Pa]` defines a lower limit below which no flow takes place (see :math:`\tau_0` at :ref:`samosatfrict`). -The dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. -:math:`\dot\gamma` is the flow velocity gradient, or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). +:math:`K` is a consistency factor :math:`[Pa \cdot s^{n}]`. For the :ref:`bingham` and the :ref:`obrienandjulien` rheologies the consistency factor corresponds to the dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` +which quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. :math:`\dot\gamma` is the flow velocity gradient, +or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). The flow index :math:`n` describes the rheological behaviour of the mixture as :cite:`KaZhHaHe2025`: - :math:`n = 0` a rate-independent solid-like behaviour, @@ -731,6 +734,8 @@ in following substitution :cite:`Iv1997, DeIv2001, GiFlSaHe2020`: where :math:`\overline{u}` is the depth-averaged flow velocity and :math:`\overline{h}` is the flow thickness. +.. _obrienandjulien: + O´Brien and Julien +++++++++++++++++++ The quadratic rheological model proposed by O´Brien and Julien :cite:`ObJu1985` accounts for various shear stress components, including @@ -761,17 +766,21 @@ of sediment :math:`[kg \cdot m^{-3}]`, :math:`d_s` is the sediment size :math:`[ The quadratic rheological model is particularly suitable for simulating debris flows with high sediment concentrations, in which energy dissipation and resistance due to turbulence and particle collisions are significant. +.. _herschelandbulkley: + Herschel and Bulkley +++++++++++++++++++++++++ The model by Herschel and Bulkley :cite:`HeBu1926` is expressed by an empirical power-law equation: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot{\gamma}^n + \tau^{(b)} = \tau_y + K \cdot \dot{\gamma}^n :label: herschelAndBulkley where the flow index :math:`n` describes the rheological behaviour of the mixture (see above). The factor in front of the shear rate -was originally introduced as a consistency factor :math:`K`. In :py:mod:`com1DFA` :math:`K` equals to the bulk dynamic viscosity :math:`\eta_m`. -This rheology applies to fine-grained soil-water mixtures that exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. +was originally introduced as a consistency factor :math:`K`.This rheology applies to fine-grained soil-water mixtures that +exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. + +.. _bingham: Bingham ++++++++++++++++ From ea63c91185c7c73dba682dbb146a013df92cf7ba Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Fri, 13 Mar 2026 09:58:43 +0100 Subject: [PATCH 6/8] Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py --- avaframe/com1DFA/DFAfunctionsCython.pyx | 9 +++----- avaframe/com1DFA/checkCfg.py | 3 +-- avaframe/com1DFA/com1DFA.py | 3 +-- avaframe/com1DFA/com1DFACfg.ini | 8 +++----- avaframe/tests/test_DFAfunctionsCython.py | 11 +++------- docs/theoryCom1DFA.rst | 25 +++++++++++++++-------- 6 files changed, 28 insertions(+), 31 deletions(-) diff --git a/avaframe/com1DFA/DFAfunctionsCython.pyx b/avaframe/com1DFA/DFAfunctionsCython.pyx index 72b859b23..3ef93f80b 100644 --- a/avaframe/com1DFA/DFAfunctionsCython.pyx +++ b/avaframe/com1DFA/DFAfunctionsCython.pyx @@ -95,10 +95,9 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double alpha2TauyObrienAndJulien = cfg.getfloat('alpha2TauyObrienAndJulien') cdef double beta2TauyObrienAndJulien = cfg.getfloat('beta2TauyObrienAndJulien') cdef double alphaObrienAndJulien = cfg.getfloat('alphaObrienAndJulien') - cdef double alpha1EtaHerschelAndBulkley = cfg.getfloat('alpha1EtaHerschelAndBulkley') - cdef double beta1EtaHerschelAndBulkley = cfg.getfloat('beta1EtaHerschelAndBulkley') cdef double alpha2TauyHerschelAndBulkley = cfg.getfloat('alpha2TauyHerschelAndBulkley') cdef double beta2TauyHerschelAndBulkley = cfg.getfloat('beta2TauyHerschelAndBulkley') + cdef double kHerschelAndBulkley = cfg.getfloat('kHerschelAndBulkley') cdef double nHerschelAndBulkley = cfg.getfloat('nHerschelAndBulkley') cdef double alpha1EtaBingham = cfg.getfloat('alpha1EtaBingham') cdef double beta1EtaBingham = cfg.getfloat('beta1EtaBingham') @@ -170,7 +169,7 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double gravAccNorm, accNormCurv, effAccNorm, gravAccTangX, gravAccTangY, gravAccTangZ, forceBotTang, sigmaB, tau cdef double muVoellmyRaster, xsiVoellmyRaster cdef double shearRate, etaObrienAndJulien, tauyObrienAndJulien, lmObrienAndJulien - cdef double lambdaBagnold, cObrienAndJulien, etaHerschelAndBulkley, tauyHerschelAndBulkley, etaBingham, tauyBingham + cdef double lambdaBagnold, cObrienAndJulien, tauyHerschelAndBulkley, etaBingham, tauyBingham # variables for interpolation cdef int Lx0, Ly0, LxEnd0, LyEnd0, iCell, iCellEnd cdef double w[4] @@ -340,12 +339,10 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType tau = tauyObrienAndJulien + etaObrienAndJulien * shearRate + cObrienAndJulien * (shearRate * shearRate) elif frictType == 11: ## Herschel and Bulkley - # viscosity - etaHerschelAndBulkley = alpha1EtaHerschelAndBulkley * math.exp(beta1EtaHerschelAndBulkley * cvSediment) # yield shear stress tauyHerschelAndBulkley = alpha2TauyHerschelAndBulkley * math.exp(beta2TauyHerschelAndBulkley * cvSediment) # shear stress - tau = tauyHerschelAndBulkley + etaHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) + tau = tauyHerschelAndBulkley + kHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) elif frictType == 12: ## Bingham # viscosity diff --git a/avaframe/com1DFA/checkCfg.py b/avaframe/com1DFA/checkCfg.py index ef8289a43..83c76f780 100644 --- a/avaframe/com1DFA/checkCfg.py +++ b/avaframe/com1DFA/checkCfg.py @@ -135,10 +135,9 @@ def checkCfgFrictionModel(cfg, inputSimFiles, relVolume=""): "alpha2TauyObrienAndJulien", "beta2TauyObrienAndJulien", "alphaObrienAndJulien", - "alpha1EtaHerschelAndBulkley", - "beta1EtaHerschelAndBulkley", "alpha2TauyHerschelAndBulkley", "beta2TauyHerschelAndBulkley", + "kHerschelAndBulkley", "nHerschelAndBulkley", "alpha1EtaBingham", "beta1EtaBingham", diff --git a/avaframe/com1DFA/com1DFA.py b/avaframe/com1DFA/com1DFA.py index ef2993f9a..20cad75a2 100644 --- a/avaframe/com1DFA/com1DFA.py +++ b/avaframe/com1DFA/com1DFA.py @@ -932,8 +932,7 @@ def createReportDict(avaDir, logName, relName, inputSimLines, cfg, reportAreaInf "type": "columns", "model": "Herschel and Bulkley", "n": cfgGen["nHerschelAndBulkley"], - "alpha1": cfgGen["alpha1EtaHerschelAndBulkley"], - "beta1": cfgGen["beta1EtaHerschelAndBulkley"], + "K": cfgGen["kHerschelAndBulkley"], "alpha2": cfgGen["alpha2TauyHerschelAndBulkley"], "beta2": cfgGen["beta2TauyHerschelAndBulkley"], } diff --git a/avaframe/com1DFA/com1DFACfg.ini b/avaframe/com1DFA/com1DFACfg.ini index 1a47466e1..dae998349 100644 --- a/avaframe/com1DFA/com1DFACfg.ini +++ b/avaframe/com1DFA/com1DFACfg.ini @@ -401,16 +401,14 @@ beta2TauyObrienAndJulien = 13.11 # parameter, default value Takahashi (1978) alphaObrienAndJulien = 0.01 #++++++++++++Herschel and Bulkley -# viscosity parameter -## TODO adjusting, default value O´Brien and Julien (1985) -alpha1EtaHerschelAndBulkley = 0.650 - ## TODO adjusting, default value O´Brien and Julien (1985) -beta1EtaHerschelAndBulkley = 16.81 # yield shear stress parameter ## TODO adjusting, default value O´Brien and Julien (1985) alpha2TauyHerschelAndBulkley = 0.00886 ## TODO adjusting, default value O´Brien and Julien (1985) beta2TauyHerschelAndBulkley = 13.11 +# consistency factor K +## TODO adjusting +kHerschelAndBulkley = 50 # flow exponent ## TODO adjusting nHerschelAndBulkley = 1.5 diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 816ce665b..731dec875 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -839,10 +839,9 @@ def test_computeForceC(): "mucoulomb": "0.155", "mucoulombminshear": "0.155", "tau0coulombminshear": "70", - "alpha1EtaHerschelAndBulkley": "0.650", - "beta1EtaHerschelAndBulkley": "16.81", "alpha2TauyHerschelAndBulkley": "0.00886", "beta2TauyHerschelAndBulkley": "13.11", + "kHerschelAndBulkley": "50", "nHerschelAndBulkley": "1.5", "alpha1EtaBingham": "0.650", "beta1EtaBingham": "16.81", @@ -985,19 +984,15 @@ def test_computeForceC(): ## Herschel and Bulkley # read input parameters - test_alpha1EtaHerschelAndBulkley = float(cfg["GENERAL"]["alpha1EtaHerschelAndBulkley"]) - test_beta1EtaHerschelAndBulkley = float(cfg["GENERAL"]["beta1EtaHerschelAndBulkley"]) test_alpha2TauyHerschelAndBulkley = float(cfg["GENERAL"]["alpha2TauyHerschelAndBulkley"]) test_beta2TauyHerschelAndBulkley = float(cfg["GENERAL"]["beta2TauyHerschelAndBulkley"]) + test_kHerschelAndBulkley = float(cfg["GENERAL"]["kHerschelAndBulkley"]) test_nHerschelAndBulkley = float(cfg["GENERAL"]["nHerschelAndBulkley"]) # compute tau - test_etaHerschelandBulkley = test_alpha1EtaHerschelAndBulkley * math.exp( - test_beta1EtaHerschelAndBulkley * test_cvSediment - ) test_tauyHerschelandBulkley = test_alpha2TauyHerschelAndBulkley * math.exp( test_beta2TauyHerschelAndBulkley * test_cvSediment ) - test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_etaHerschelandBulkley * math.pow( + test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_kHerschelAndBulkley * math.pow( test_shearRate[0], test_nHerschelAndBulkley ) test_tauArray[1] = test_tauHerschelAndBulkley diff --git a/docs/theoryCom1DFA.rst b/docs/theoryCom1DFA.rst index da3fce612..ceca2977d 100644 --- a/docs/theoryCom1DFA.rst +++ b/docs/theoryCom1DFA.rst @@ -480,10 +480,11 @@ With &\mathbf{x} \qquad &\text{position vector}\end{aligned} For solid-fluid mixtures like debris flows, where the properties of the fluid phase are dominating the flow process, -the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). +the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` (resp. consistency factor :math:`K` for :ref:`herschelandbulkley`) +and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). .. math:: - \tau^{(b)}_i = f(\eta_m,\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) + \tau^{(b)}_i = f(\eta_m,(K),\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) :label: rheological model With @@ -491,6 +492,7 @@ With .. math:: \begin{aligned} &\eta_m \qquad &\text{dynamic viscosity of the solid-fluid mixture}\\ + &K \qquad &\text{consistency factor of the solid-fluid mixture}\\ &\dot\gamma \qquad &\text{shear rate}\\ &\overline{h} \qquad &\text{average flow thickness}\\ &\rho_m \qquad &\text{density of the solid-fluid mixture}\\ @@ -662,12 +664,13 @@ and non-Newtonian fluids, in which a yield shear stress has to be exceeded or sh The rheological models are incorporated into a single, general form, which can be expressed as follows: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot\gamma^n + C \cdot \dot\gamma^2 + \tau^{(b)} = \tau_y + K \cdot \dot\gamma^n + C \cdot \dot\gamma^2 :label: rheology-general The yield shear stress :math:`\tau_y` :math:`[Pa]` defines a lower limit below which no flow takes place (see :math:`\tau_0` at :ref:`samosatfrict`). -The dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. -:math:`\dot\gamma` is the flow velocity gradient, or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). +:math:`K` is a consistency factor :math:`[Pa \cdot s^{n}]`. For the :ref:`bingham` and the :ref:`obrienandjulien` rheologies the consistency factor corresponds to the dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` +which quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. :math:`\dot\gamma` is the flow velocity gradient, +or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). The flow index :math:`n` describes the rheological behaviour of the mixture as :cite:`KaZhHaHe2025`: - :math:`n = 0` a rate-independent solid-like behaviour, @@ -731,6 +734,8 @@ in following substitution :cite:`Iv1997, DeIv2001, GiFlSaHe2020`: where :math:`\overline{u}` is the depth-averaged flow velocity and :math:`\overline{h}` is the flow thickness. +.. _obrienandjulien: + O´Brien and Julien +++++++++++++++++++ The quadratic rheological model proposed by O´Brien and Julien :cite:`ObJu1985` accounts for various shear stress components, including @@ -761,17 +766,21 @@ of sediment :math:`[kg \cdot m^{-3}]`, :math:`d_s` is the sediment size :math:`[ The quadratic rheological model is particularly suitable for simulating debris flows with high sediment concentrations, in which energy dissipation and resistance due to turbulence and particle collisions are significant. +.. _herschelandbulkley: + Herschel and Bulkley +++++++++++++++++++++++++ The model by Herschel and Bulkley :cite:`HeBu1926` is expressed by an empirical power-law equation: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot{\gamma}^n + \tau^{(b)} = \tau_y + K \cdot \dot{\gamma}^n :label: herschelAndBulkley where the flow index :math:`n` describes the rheological behaviour of the mixture (see above). The factor in front of the shear rate -was originally introduced as a consistency factor :math:`K`. In :py:mod:`com1DFA` :math:`K` equals to the bulk dynamic viscosity :math:`\eta_m`. -This rheology applies to fine-grained soil-water mixtures that exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. +was originally introduced as a consistency factor :math:`K`.This rheology applies to fine-grained soil-water mixtures that +exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. + +.. _bingham: Bingham ++++++++++++++++ From 483c5057f6b8233dd07c9fd84a1ef9db5ef84ab7 Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Fri, 13 Mar 2026 09:58:43 +0100 Subject: [PATCH 7/8] Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py --- avaframe/com1DFA/DFAfunctionsCython.pyx | 9 +++----- avaframe/com1DFA/checkCfg.py | 3 +-- avaframe/com1DFA/com1DFA.py | 3 +-- avaframe/com1DFA/com1DFACfg.ini | 8 +++----- avaframe/tests/test_DFAfunctionsCython.py | 11 +++------- docs/theoryCom1DFA.rst | 25 +++++++++++++++-------- 6 files changed, 28 insertions(+), 31 deletions(-) diff --git a/avaframe/com1DFA/DFAfunctionsCython.pyx b/avaframe/com1DFA/DFAfunctionsCython.pyx index 72b859b23..3ef93f80b 100644 --- a/avaframe/com1DFA/DFAfunctionsCython.pyx +++ b/avaframe/com1DFA/DFAfunctionsCython.pyx @@ -95,10 +95,9 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double alpha2TauyObrienAndJulien = cfg.getfloat('alpha2TauyObrienAndJulien') cdef double beta2TauyObrienAndJulien = cfg.getfloat('beta2TauyObrienAndJulien') cdef double alphaObrienAndJulien = cfg.getfloat('alphaObrienAndJulien') - cdef double alpha1EtaHerschelAndBulkley = cfg.getfloat('alpha1EtaHerschelAndBulkley') - cdef double beta1EtaHerschelAndBulkley = cfg.getfloat('beta1EtaHerschelAndBulkley') cdef double alpha2TauyHerschelAndBulkley = cfg.getfloat('alpha2TauyHerschelAndBulkley') cdef double beta2TauyHerschelAndBulkley = cfg.getfloat('beta2TauyHerschelAndBulkley') + cdef double kHerschelAndBulkley = cfg.getfloat('kHerschelAndBulkley') cdef double nHerschelAndBulkley = cfg.getfloat('nHerschelAndBulkley') cdef double alpha1EtaBingham = cfg.getfloat('alpha1EtaBingham') cdef double beta1EtaBingham = cfg.getfloat('beta1EtaBingham') @@ -170,7 +169,7 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double gravAccNorm, accNormCurv, effAccNorm, gravAccTangX, gravAccTangY, gravAccTangZ, forceBotTang, sigmaB, tau cdef double muVoellmyRaster, xsiVoellmyRaster cdef double shearRate, etaObrienAndJulien, tauyObrienAndJulien, lmObrienAndJulien - cdef double lambdaBagnold, cObrienAndJulien, etaHerschelAndBulkley, tauyHerschelAndBulkley, etaBingham, tauyBingham + cdef double lambdaBagnold, cObrienAndJulien, tauyHerschelAndBulkley, etaBingham, tauyBingham # variables for interpolation cdef int Lx0, Ly0, LxEnd0, LyEnd0, iCell, iCellEnd cdef double w[4] @@ -340,12 +339,10 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType tau = tauyObrienAndJulien + etaObrienAndJulien * shearRate + cObrienAndJulien * (shearRate * shearRate) elif frictType == 11: ## Herschel and Bulkley - # viscosity - etaHerschelAndBulkley = alpha1EtaHerschelAndBulkley * math.exp(beta1EtaHerschelAndBulkley * cvSediment) # yield shear stress tauyHerschelAndBulkley = alpha2TauyHerschelAndBulkley * math.exp(beta2TauyHerschelAndBulkley * cvSediment) # shear stress - tau = tauyHerschelAndBulkley + etaHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) + tau = tauyHerschelAndBulkley + kHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) elif frictType == 12: ## Bingham # viscosity diff --git a/avaframe/com1DFA/checkCfg.py b/avaframe/com1DFA/checkCfg.py index ef8289a43..83c76f780 100644 --- a/avaframe/com1DFA/checkCfg.py +++ b/avaframe/com1DFA/checkCfg.py @@ -135,10 +135,9 @@ def checkCfgFrictionModel(cfg, inputSimFiles, relVolume=""): "alpha2TauyObrienAndJulien", "beta2TauyObrienAndJulien", "alphaObrienAndJulien", - "alpha1EtaHerschelAndBulkley", - "beta1EtaHerschelAndBulkley", "alpha2TauyHerschelAndBulkley", "beta2TauyHerschelAndBulkley", + "kHerschelAndBulkley", "nHerschelAndBulkley", "alpha1EtaBingham", "beta1EtaBingham", diff --git a/avaframe/com1DFA/com1DFA.py b/avaframe/com1DFA/com1DFA.py index ef2993f9a..20cad75a2 100644 --- a/avaframe/com1DFA/com1DFA.py +++ b/avaframe/com1DFA/com1DFA.py @@ -932,8 +932,7 @@ def createReportDict(avaDir, logName, relName, inputSimLines, cfg, reportAreaInf "type": "columns", "model": "Herschel and Bulkley", "n": cfgGen["nHerschelAndBulkley"], - "alpha1": cfgGen["alpha1EtaHerschelAndBulkley"], - "beta1": cfgGen["beta1EtaHerschelAndBulkley"], + "K": cfgGen["kHerschelAndBulkley"], "alpha2": cfgGen["alpha2TauyHerschelAndBulkley"], "beta2": cfgGen["beta2TauyHerschelAndBulkley"], } diff --git a/avaframe/com1DFA/com1DFACfg.ini b/avaframe/com1DFA/com1DFACfg.ini index 1a47466e1..dae998349 100644 --- a/avaframe/com1DFA/com1DFACfg.ini +++ b/avaframe/com1DFA/com1DFACfg.ini @@ -401,16 +401,14 @@ beta2TauyObrienAndJulien = 13.11 # parameter, default value Takahashi (1978) alphaObrienAndJulien = 0.01 #++++++++++++Herschel and Bulkley -# viscosity parameter -## TODO adjusting, default value O´Brien and Julien (1985) -alpha1EtaHerschelAndBulkley = 0.650 - ## TODO adjusting, default value O´Brien and Julien (1985) -beta1EtaHerschelAndBulkley = 16.81 # yield shear stress parameter ## TODO adjusting, default value O´Brien and Julien (1985) alpha2TauyHerschelAndBulkley = 0.00886 ## TODO adjusting, default value O´Brien and Julien (1985) beta2TauyHerschelAndBulkley = 13.11 +# consistency factor K +## TODO adjusting +kHerschelAndBulkley = 50 # flow exponent ## TODO adjusting nHerschelAndBulkley = 1.5 diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 816ce665b..731dec875 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -839,10 +839,9 @@ def test_computeForceC(): "mucoulomb": "0.155", "mucoulombminshear": "0.155", "tau0coulombminshear": "70", - "alpha1EtaHerschelAndBulkley": "0.650", - "beta1EtaHerschelAndBulkley": "16.81", "alpha2TauyHerschelAndBulkley": "0.00886", "beta2TauyHerschelAndBulkley": "13.11", + "kHerschelAndBulkley": "50", "nHerschelAndBulkley": "1.5", "alpha1EtaBingham": "0.650", "beta1EtaBingham": "16.81", @@ -985,19 +984,15 @@ def test_computeForceC(): ## Herschel and Bulkley # read input parameters - test_alpha1EtaHerschelAndBulkley = float(cfg["GENERAL"]["alpha1EtaHerschelAndBulkley"]) - test_beta1EtaHerschelAndBulkley = float(cfg["GENERAL"]["beta1EtaHerschelAndBulkley"]) test_alpha2TauyHerschelAndBulkley = float(cfg["GENERAL"]["alpha2TauyHerschelAndBulkley"]) test_beta2TauyHerschelAndBulkley = float(cfg["GENERAL"]["beta2TauyHerschelAndBulkley"]) + test_kHerschelAndBulkley = float(cfg["GENERAL"]["kHerschelAndBulkley"]) test_nHerschelAndBulkley = float(cfg["GENERAL"]["nHerschelAndBulkley"]) # compute tau - test_etaHerschelandBulkley = test_alpha1EtaHerschelAndBulkley * math.exp( - test_beta1EtaHerschelAndBulkley * test_cvSediment - ) test_tauyHerschelandBulkley = test_alpha2TauyHerschelAndBulkley * math.exp( test_beta2TauyHerschelAndBulkley * test_cvSediment ) - test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_etaHerschelandBulkley * math.pow( + test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_kHerschelAndBulkley * math.pow( test_shearRate[0], test_nHerschelAndBulkley ) test_tauArray[1] = test_tauHerschelAndBulkley diff --git a/docs/theoryCom1DFA.rst b/docs/theoryCom1DFA.rst index da3fce612..ceca2977d 100644 --- a/docs/theoryCom1DFA.rst +++ b/docs/theoryCom1DFA.rst @@ -480,10 +480,11 @@ With &\mathbf{x} \qquad &\text{position vector}\end{aligned} For solid-fluid mixtures like debris flows, where the properties of the fluid phase are dominating the flow process, -the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). +the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` (resp. consistency factor :math:`K` for :ref:`herschelandbulkley`) +and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). .. math:: - \tau^{(b)}_i = f(\eta_m,\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) + \tau^{(b)}_i = f(\eta_m,(K),\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) :label: rheological model With @@ -491,6 +492,7 @@ With .. math:: \begin{aligned} &\eta_m \qquad &\text{dynamic viscosity of the solid-fluid mixture}\\ + &K \qquad &\text{consistency factor of the solid-fluid mixture}\\ &\dot\gamma \qquad &\text{shear rate}\\ &\overline{h} \qquad &\text{average flow thickness}\\ &\rho_m \qquad &\text{density of the solid-fluid mixture}\\ @@ -662,12 +664,13 @@ and non-Newtonian fluids, in which a yield shear stress has to be exceeded or sh The rheological models are incorporated into a single, general form, which can be expressed as follows: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot\gamma^n + C \cdot \dot\gamma^2 + \tau^{(b)} = \tau_y + K \cdot \dot\gamma^n + C \cdot \dot\gamma^2 :label: rheology-general The yield shear stress :math:`\tau_y` :math:`[Pa]` defines a lower limit below which no flow takes place (see :math:`\tau_0` at :ref:`samosatfrict`). -The dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. -:math:`\dot\gamma` is the flow velocity gradient, or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). +:math:`K` is a consistency factor :math:`[Pa \cdot s^{n}]`. For the :ref:`bingham` and the :ref:`obrienandjulien` rheologies the consistency factor corresponds to the dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` +which quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. :math:`\dot\gamma` is the flow velocity gradient, +or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). The flow index :math:`n` describes the rheological behaviour of the mixture as :cite:`KaZhHaHe2025`: - :math:`n = 0` a rate-independent solid-like behaviour, @@ -731,6 +734,8 @@ in following substitution :cite:`Iv1997, DeIv2001, GiFlSaHe2020`: where :math:`\overline{u}` is the depth-averaged flow velocity and :math:`\overline{h}` is the flow thickness. +.. _obrienandjulien: + O´Brien and Julien +++++++++++++++++++ The quadratic rheological model proposed by O´Brien and Julien :cite:`ObJu1985` accounts for various shear stress components, including @@ -761,17 +766,21 @@ of sediment :math:`[kg \cdot m^{-3}]`, :math:`d_s` is the sediment size :math:`[ The quadratic rheological model is particularly suitable for simulating debris flows with high sediment concentrations, in which energy dissipation and resistance due to turbulence and particle collisions are significant. +.. _herschelandbulkley: + Herschel and Bulkley +++++++++++++++++++++++++ The model by Herschel and Bulkley :cite:`HeBu1926` is expressed by an empirical power-law equation: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot{\gamma}^n + \tau^{(b)} = \tau_y + K \cdot \dot{\gamma}^n :label: herschelAndBulkley where the flow index :math:`n` describes the rheological behaviour of the mixture (see above). The factor in front of the shear rate -was originally introduced as a consistency factor :math:`K`. In :py:mod:`com1DFA` :math:`K` equals to the bulk dynamic viscosity :math:`\eta_m`. -This rheology applies to fine-grained soil-water mixtures that exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. +was originally introduced as a consistency factor :math:`K`.This rheology applies to fine-grained soil-water mixtures that +exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. + +.. _bingham: Bingham ++++++++++++++++ From 85bb4b40d558d78dce289bc56bdfe2ee047d7626 Mon Sep 17 00:00:00 2001 From: JuLa96 <206707905+JuLa96@users.noreply.github.com> Date: Fri, 13 Mar 2026 09:58:43 +0100 Subject: [PATCH 8/8] Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py Replace viscosity by consistency factor for Herschel-Bulkley rheology Replace viscosity by consistency factor in checkCfg.py Fix input consistency factor in avaframe/tests/test_DFAfunctionsCython.py --- avaframe/com1DFA/DFAfunctionsCython.pyx | 9 +++----- avaframe/com1DFA/checkCfg.py | 3 +-- avaframe/com1DFA/com1DFA.py | 3 +-- avaframe/com1DFA/com1DFACfg.ini | 8 +++----- avaframe/tests/test_DFAfunctionsCython.py | 11 +++------- docs/theoryCom1DFA.rst | 25 +++++++++++++++-------- 6 files changed, 28 insertions(+), 31 deletions(-) diff --git a/avaframe/com1DFA/DFAfunctionsCython.pyx b/avaframe/com1DFA/DFAfunctionsCython.pyx index 72b859b23..3ef93f80b 100644 --- a/avaframe/com1DFA/DFAfunctionsCython.pyx +++ b/avaframe/com1DFA/DFAfunctionsCython.pyx @@ -95,10 +95,9 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double alpha2TauyObrienAndJulien = cfg.getfloat('alpha2TauyObrienAndJulien') cdef double beta2TauyObrienAndJulien = cfg.getfloat('beta2TauyObrienAndJulien') cdef double alphaObrienAndJulien = cfg.getfloat('alphaObrienAndJulien') - cdef double alpha1EtaHerschelAndBulkley = cfg.getfloat('alpha1EtaHerschelAndBulkley') - cdef double beta1EtaHerschelAndBulkley = cfg.getfloat('beta1EtaHerschelAndBulkley') cdef double alpha2TauyHerschelAndBulkley = cfg.getfloat('alpha2TauyHerschelAndBulkley') cdef double beta2TauyHerschelAndBulkley = cfg.getfloat('beta2TauyHerschelAndBulkley') + cdef double kHerschelAndBulkley = cfg.getfloat('kHerschelAndBulkley') cdef double nHerschelAndBulkley = cfg.getfloat('nHerschelAndBulkley') cdef double alpha1EtaBingham = cfg.getfloat('alpha1EtaBingham') cdef double beta1EtaBingham = cfg.getfloat('beta1EtaBingham') @@ -170,7 +169,7 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType cdef double gravAccNorm, accNormCurv, effAccNorm, gravAccTangX, gravAccTangY, gravAccTangZ, forceBotTang, sigmaB, tau cdef double muVoellmyRaster, xsiVoellmyRaster cdef double shearRate, etaObrienAndJulien, tauyObrienAndJulien, lmObrienAndJulien - cdef double lambdaBagnold, cObrienAndJulien, etaHerschelAndBulkley, tauyHerschelAndBulkley, etaBingham, tauyBingham + cdef double lambdaBagnold, cObrienAndJulien, tauyHerschelAndBulkley, etaBingham, tauyBingham # variables for interpolation cdef int Lx0, Ly0, LxEnd0, LyEnd0, iCell, iCellEnd cdef double w[4] @@ -340,12 +339,10 @@ def computeForceC(cfg, particles, fields, dem, int frictType, int resistanceType tau = tauyObrienAndJulien + etaObrienAndJulien * shearRate + cObrienAndJulien * (shearRate * shearRate) elif frictType == 11: ## Herschel and Bulkley - # viscosity - etaHerschelAndBulkley = alpha1EtaHerschelAndBulkley * math.exp(beta1EtaHerschelAndBulkley * cvSediment) # yield shear stress tauyHerschelAndBulkley = alpha2TauyHerschelAndBulkley * math.exp(beta2TauyHerschelAndBulkley * cvSediment) # shear stress - tau = tauyHerschelAndBulkley + etaHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) + tau = tauyHerschelAndBulkley + kHerschelAndBulkley * math.pow(shearRate, nHerschelAndBulkley) elif frictType == 12: ## Bingham # viscosity diff --git a/avaframe/com1DFA/checkCfg.py b/avaframe/com1DFA/checkCfg.py index ef8289a43..83c76f780 100644 --- a/avaframe/com1DFA/checkCfg.py +++ b/avaframe/com1DFA/checkCfg.py @@ -135,10 +135,9 @@ def checkCfgFrictionModel(cfg, inputSimFiles, relVolume=""): "alpha2TauyObrienAndJulien", "beta2TauyObrienAndJulien", "alphaObrienAndJulien", - "alpha1EtaHerschelAndBulkley", - "beta1EtaHerschelAndBulkley", "alpha2TauyHerschelAndBulkley", "beta2TauyHerschelAndBulkley", + "kHerschelAndBulkley", "nHerschelAndBulkley", "alpha1EtaBingham", "beta1EtaBingham", diff --git a/avaframe/com1DFA/com1DFA.py b/avaframe/com1DFA/com1DFA.py index ef2993f9a..20cad75a2 100644 --- a/avaframe/com1DFA/com1DFA.py +++ b/avaframe/com1DFA/com1DFA.py @@ -932,8 +932,7 @@ def createReportDict(avaDir, logName, relName, inputSimLines, cfg, reportAreaInf "type": "columns", "model": "Herschel and Bulkley", "n": cfgGen["nHerschelAndBulkley"], - "alpha1": cfgGen["alpha1EtaHerschelAndBulkley"], - "beta1": cfgGen["beta1EtaHerschelAndBulkley"], + "K": cfgGen["kHerschelAndBulkley"], "alpha2": cfgGen["alpha2TauyHerschelAndBulkley"], "beta2": cfgGen["beta2TauyHerschelAndBulkley"], } diff --git a/avaframe/com1DFA/com1DFACfg.ini b/avaframe/com1DFA/com1DFACfg.ini index 1a47466e1..dae998349 100644 --- a/avaframe/com1DFA/com1DFACfg.ini +++ b/avaframe/com1DFA/com1DFACfg.ini @@ -401,16 +401,14 @@ beta2TauyObrienAndJulien = 13.11 # parameter, default value Takahashi (1978) alphaObrienAndJulien = 0.01 #++++++++++++Herschel and Bulkley -# viscosity parameter -## TODO adjusting, default value O´Brien and Julien (1985) -alpha1EtaHerschelAndBulkley = 0.650 - ## TODO adjusting, default value O´Brien and Julien (1985) -beta1EtaHerschelAndBulkley = 16.81 # yield shear stress parameter ## TODO adjusting, default value O´Brien and Julien (1985) alpha2TauyHerschelAndBulkley = 0.00886 ## TODO adjusting, default value O´Brien and Julien (1985) beta2TauyHerschelAndBulkley = 13.11 +# consistency factor K +## TODO adjusting +kHerschelAndBulkley = 50 # flow exponent ## TODO adjusting nHerschelAndBulkley = 1.5 diff --git a/avaframe/tests/test_DFAfunctionsCython.py b/avaframe/tests/test_DFAfunctionsCython.py index 816ce665b..731dec875 100644 --- a/avaframe/tests/test_DFAfunctionsCython.py +++ b/avaframe/tests/test_DFAfunctionsCython.py @@ -839,10 +839,9 @@ def test_computeForceC(): "mucoulomb": "0.155", "mucoulombminshear": "0.155", "tau0coulombminshear": "70", - "alpha1EtaHerschelAndBulkley": "0.650", - "beta1EtaHerschelAndBulkley": "16.81", "alpha2TauyHerschelAndBulkley": "0.00886", "beta2TauyHerschelAndBulkley": "13.11", + "kHerschelAndBulkley": "50", "nHerschelAndBulkley": "1.5", "alpha1EtaBingham": "0.650", "beta1EtaBingham": "16.81", @@ -985,19 +984,15 @@ def test_computeForceC(): ## Herschel and Bulkley # read input parameters - test_alpha1EtaHerschelAndBulkley = float(cfg["GENERAL"]["alpha1EtaHerschelAndBulkley"]) - test_beta1EtaHerschelAndBulkley = float(cfg["GENERAL"]["beta1EtaHerschelAndBulkley"]) test_alpha2TauyHerschelAndBulkley = float(cfg["GENERAL"]["alpha2TauyHerschelAndBulkley"]) test_beta2TauyHerschelAndBulkley = float(cfg["GENERAL"]["beta2TauyHerschelAndBulkley"]) + test_kHerschelAndBulkley = float(cfg["GENERAL"]["kHerschelAndBulkley"]) test_nHerschelAndBulkley = float(cfg["GENERAL"]["nHerschelAndBulkley"]) # compute tau - test_etaHerschelandBulkley = test_alpha1EtaHerschelAndBulkley * math.exp( - test_beta1EtaHerschelAndBulkley * test_cvSediment - ) test_tauyHerschelandBulkley = test_alpha2TauyHerschelAndBulkley * math.exp( test_beta2TauyHerschelAndBulkley * test_cvSediment ) - test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_etaHerschelandBulkley * math.pow( + test_tauHerschelAndBulkley = test_tauyHerschelandBulkley + test_kHerschelAndBulkley * math.pow( test_shearRate[0], test_nHerschelAndBulkley ) test_tauArray[1] = test_tauHerschelAndBulkley diff --git a/docs/theoryCom1DFA.rst b/docs/theoryCom1DFA.rst index da3fce612..ceca2977d 100644 --- a/docs/theoryCom1DFA.rst +++ b/docs/theoryCom1DFA.rst @@ -480,10 +480,11 @@ With &\mathbf{x} \qquad &\text{position vector}\end{aligned} For solid-fluid mixtures like debris flows, where the properties of the fluid phase are dominating the flow process, -the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). +the shear stress tensor is, among other things, a function of dynamic viscosity :math:`\eta_m` (resp. consistency factor :math:`K` for :ref:`herschelandbulkley`) +and shear rate :math:`\dot\gamma` (:ref:`theoryCom1DFA:Rheological models`). .. math:: - \tau^{(b)}_i = f(\eta_m,\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) + \tau^{(b)}_i = f(\eta_m,(K),\dot\gamma,\overline{h},\rho_m,t,\mathbf{x}) :label: rheological model With @@ -491,6 +492,7 @@ With .. math:: \begin{aligned} &\eta_m \qquad &\text{dynamic viscosity of the solid-fluid mixture}\\ + &K \qquad &\text{consistency factor of the solid-fluid mixture}\\ &\dot\gamma \qquad &\text{shear rate}\\ &\overline{h} \qquad &\text{average flow thickness}\\ &\rho_m \qquad &\text{density of the solid-fluid mixture}\\ @@ -662,12 +664,13 @@ and non-Newtonian fluids, in which a yield shear stress has to be exceeded or sh The rheological models are incorporated into a single, general form, which can be expressed as follows: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot\gamma^n + C \cdot \dot\gamma^2 + \tau^{(b)} = \tau_y + K \cdot \dot\gamma^n + C \cdot \dot\gamma^2 :label: rheology-general The yield shear stress :math:`\tau_y` :math:`[Pa]` defines a lower limit below which no flow takes place (see :math:`\tau_0` at :ref:`samosatfrict`). -The dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. -:math:`\dot\gamma` is the flow velocity gradient, or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). +:math:`K` is a consistency factor :math:`[Pa \cdot s^{n}]`. For the :ref:`bingham` and the :ref:`obrienandjulien` rheologies the consistency factor corresponds to the dynamic viscosity :math:`\eta_m` :math:`[Pa \cdot s]` +which quantifies the internal frictional force between two neighbouring layers of the mixture in relative motion. :math:`\dot\gamma` is the flow velocity gradient, +or shear rate, :math:`\frac{du}{dz}` :math:`[s^{-1}]` along the axis normal to the slope (flow thickness). The flow index :math:`n` describes the rheological behaviour of the mixture as :cite:`KaZhHaHe2025`: - :math:`n = 0` a rate-independent solid-like behaviour, @@ -731,6 +734,8 @@ in following substitution :cite:`Iv1997, DeIv2001, GiFlSaHe2020`: where :math:`\overline{u}` is the depth-averaged flow velocity and :math:`\overline{h}` is the flow thickness. +.. _obrienandjulien: + O´Brien and Julien +++++++++++++++++++ The quadratic rheological model proposed by O´Brien and Julien :cite:`ObJu1985` accounts for various shear stress components, including @@ -761,17 +766,21 @@ of sediment :math:`[kg \cdot m^{-3}]`, :math:`d_s` is the sediment size :math:`[ The quadratic rheological model is particularly suitable for simulating debris flows with high sediment concentrations, in which energy dissipation and resistance due to turbulence and particle collisions are significant. +.. _herschelandbulkley: + Herschel and Bulkley +++++++++++++++++++++++++ The model by Herschel and Bulkley :cite:`HeBu1926` is expressed by an empirical power-law equation: .. math:: - \tau^{(b)} = \tau_y + \eta_m \cdot \dot{\gamma}^n + \tau^{(b)} = \tau_y + K \cdot \dot{\gamma}^n :label: herschelAndBulkley where the flow index :math:`n` describes the rheological behaviour of the mixture (see above). The factor in front of the shear rate -was originally introduced as a consistency factor :math:`K`. In :py:mod:`com1DFA` :math:`K` equals to the bulk dynamic viscosity :math:`\eta_m`. -This rheology applies to fine-grained soil-water mixtures that exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. +was originally introduced as a consistency factor :math:`K`.This rheology applies to fine-grained soil-water mixtures that +exhibit shear-thinning or shear-thickening behavior, respectively, with increasing shear rates. + +.. _bingham: Bingham ++++++++++++++++