Processing math: 100%

Steady state calculation

Calculation at a section

Transport equation

The steady-state advection equation is:

dCQdx=E(x,C)S

Integration using fourth-order Runge-Kutta méthod

Previous equation is rewriten as dydx=f(y,x) (where y=CQ and f(y,x)=E(x,y)S) and integrated using a fourth-order Runge-Kutta méthod(RK4) :

yj+1=yj+Δx6(k1+2k2+2k3+k4) where :

  • k1=f(yj,xj)=E(xj,yj)Sj=E(xj,Cj)Sj
  • k2=f(yj+Δx2k1,xj+12)=E(xj+12,yj+Δx2k1)Sj+12=E(xj+12,yj+Δx2k1Qj+12)Sj+12 where Sj+12=Sj+Sj+12 and
    Qj+12=Qj+Qj+12
  • k3=f(yj+Δx2k2,xj+12)=E(xj+12,yj+Δx2k2)Sj+12=E(xj+12,yj+Δx2k2Qj+12)Sj+12
  • k4=f(yj+Δxk3,xj+1)=E(xj+1,yj+Δxk3)Sj+1=E(xj+1,yj+Δxk3Qj+1)Sj+1

Concentrations in section j+1 are therefore :

Cj+1=yj+1Qj+1

Calculation at a node

Continuity equation at a node

The steady state mass balance at a node is :

sbQb>0CbsbQb+Qp>0CpQp+VCasE(t,Cnd)+Cnd[sbQb<0kbsbQb+Qp<0kpQpkinfScasvinf]=0

Newton’s method

Let f(Cind)=0 be the mass balance equation for the node on iteration i.

Newton’s method gives :

Cind=Ci1ndf(Ci1nd)f(Ci1nd)=Ci1ndbiai

bi=sbQb>0CbsbQb+Qp>0CpQp+VCasE(t,Ci1nd)
+Ci1nd[sbQb<0kbsbQb+Qp<0kpQpkinfScasvinf]

Is the node does not contain a pond, there is no source term and the equation is linear. If distribution coefficients kb, kp and kinf are heterogeneous, they are adjusted before the first and only iteration.

If there is a pond, the equation is no longer linear and several iterations are needed but distribution coefficients are not adjusted.

The general equation for ai is :

ai=VCasE(Ci1nd)+sbQb<0kbsbQb+Qp<0kpQpkinfScasvinf

where E(Ci1nd) is the derivative of the exchange.

Calculation for a node with a pond

Some terms in the previous equation are constant and do not need to be recomputed on every iteration :

b1=sbQb>0CbsbQb+Qp>0CpQp
a1=sbQb<0kbsbQb+Qp<0kpQpkinfScasvinf

Which yields :

bi=b1+VCasE(t,Ci1nd)+a1Ci1nd
et ai=VCasE(Ci1nd)+a1

Calculation for a node without a pond

For a node without a pond, the node concentration is directly calculated as a pondered mean of inflow concentrations.

Cnd=sbQb>0CbsbQb+Qp>0CpQpsbQb>0sbQb+Qp>0Qp

If the distribution coefficients are adjustable, Ka is computed before the iterations, as described in "Network computation".