ReacTran (version 1.4.3.1)

tran.cylindrical: Diffusive Transport in cylindrical (r, theta, z) and spherical (r, theta, phi) coordinates.

Description

Estimates the transport term (i.e. the rate of change of a concentration due to diffusion) in a cylindrical (r, theta, z) or spherical (r, theta, phi) coordinate system.

Usage

tran.cylindrical (C, C.r.up = NULL, C.r.down = NULL, 
                  C.theta.up = NULL, C.theta.down = NULL, 
                  C.z.up = NULL, C.z.down = NULL, 
                  flux.r.up = NULL, flux.r.down = NULL, 
                  flux.theta.up = NULL, flux.theta.down = NULL,          
                  flux.z.up = NULL, flux.z.down = NULL, 
                  cyclicBnd = NULL,
                  D.r = NULL, D.theta = D.r, D.z = D.r, 
                  r = NULL, theta = NULL, z = NULL)
tran.spherical (C, C.r.up = NULL, C.r.down = NULL, 
                C.theta.up = NULL, C.theta.down = NULL, 
                C.phi.up = NULL, C.phi.down = NULL, 
                flux.r.up = NULL, flux.r.down = NULL, 
                flux.theta.up = NULL, flux.theta.down = NULL,          
                flux.phi.up = NULL, flux.phi.down = NULL, 
                cyclicBnd = NULL,
                D.r = NULL, D.theta = D.r, D.phi = D.r, 
                r = NULL, theta = NULL, phi = NULL)

Arguments

C

concentration, expressed per unit volume, defined at the centre of each grid cell; Nr*Nteta*Nz (cylindrica) or Nr*Ntheta*Nphi (spherical coordinates) array [M/L3].

C.r.up

concentration at upstream boundary in r(x)-direction; one value [M/L3].

C.r.down

concentration at downstream boundary in r(x)-direction; one value [M/L3].

C.theta.up

concentration at upstream boundary in theta-direction; one value [M/L3].

C.theta.down

concentration at downstream boundary in theta-direction; one value [M/L3].

C.z.up

concentration at upstream boundary in z-direction (cylindrical coordinates); one value [M/L3].

C.z.down

concentration at downstream boundary in z-direction(cylindrical coordinates); one value [M/L3].

C.phi.up

concentration at upstream boundary in phi-direction (spherical coordinates); one value [M/L3].

C.phi.down

concentration at downstream boundary in phi-direction(spherical coordinates); one value [M/L3].

flux.r.up

flux across the upstream boundary in r-direction, positive = INTO model domain; one value [M/L2/T].

flux.r.down

flux across the downstream boundary in r-direction, positive = OUT of model domain; one value [M/L2/T].

flux.theta.up

flux across the upstream boundary in theta-direction, positive = INTO model domain; one value [M/L2/T].

flux.theta.down

flux across the downstream boundary in theta-direction, positive = OUT of model domain; one value [M/L2/T].

flux.z.up

flux across the upstream boundary in z-direction(cylindrical coordinates); positive = INTO model domain; one value [M/L2/T].

flux.z.down

flux across the downstream boundary in z-direction, (cylindrical coordinates); positive = OUT of model domain; one value [M/L2/T].

flux.phi.up

flux across the upstream boundary in phi-direction(spherical coordinates); positive = INTO model domain; one value [M/L2/T].

flux.phi.down

flux across the downstream boundary in phi-direction, (spherical coordinates); positive = OUT of model domain; one value [M/L2/T].

cyclicBnd

If not NULL, the direction in which a cyclic boundary is defined, i.e. cyclicBnd = 1 for the r direction, cyclicBnd = 2 for the theta direction and cyclicBnd = c(1,2) for both the r and theta direction.

D.r

diffusion coefficient in r-direction, defined on grid cell interfaces. One value or a vector of length (Nr+1), [L2/T].

D.theta

diffusion coefficient in theta-direction, defined on grid cell interfaces. One value or or a vector of length (Ntheta+1), [L2/T].

D.z

diffusion coefficient in z-direction, defined on grid cell interfaces for cylindrical coordinates. One value or a vector of length (Nz+1) [L2/T].

D.phi

diffusion coefficient in phi-direction, defined on grid cell interfaces for cylindrical coordinates. One value or a vector of length (Nphi+1) [L2/T].

r

position of adjacent cell interfaces in the r-direction. A vector of length Nr+1 [L].

theta

position of adjacent cell interfaces in the theta-direction. A vector of length Ntheta+1 [L]. Theta should be within [0,2 pi]

z

position of adjacent cell interfaces in the z-direction (cylindrical coordinates). A vector of length Nz+1 [L].

phi

position of adjacent cell interfaces in the phi-direction (spherical coordinates). A vector of length Nphi+1 [L]. Phi should be within [0,2 pi]

Value

a list containing:

dC

the rate of change of the concentration C due to transport, defined in the centre of each grid cell, a Nr*Nteta*Nz (cylindrical) or Nr*Ntheta*Nphi (spherical coordinates) array. [M/L3/T].

flux.r.up

flux across the upstream boundary in r-direction, positive = INTO model domain. A matrix of dimension Nteta*Nz (cylindrical) or Ntheta*Nphi (spherical) [M/L2/T].

flux.r.down

flux across the downstream boundary in r-direction, positive = OUT of model domain. A matrix of dimension Nteta*Nz (cylindrical) or Ntheta*Nphi (spherical) [M/L2/T].

flux.theta.up

flux across the upstream boundary in theta-direction, positive = INTO model domain. A matrix of dimension Nr*Nz (cylindrical) or or Nr*Nphi (spherical) [M/L2/T].

flux.theta.down

flux across the downstream boundary in theta-direction, positive = OUT of model domain. A matrix of dimension Nr*Nz (cylindrical) or Nr*Nphi (spherical) [M/L2/T].

flux.z.up

flux across the upstream boundary in z-direction, for cylindrical coordinates; positive = OUT of model domain. A matrix of dimension Nr*Nteta [M/L2/T].

flux.z.down

flux across the downstream boundary in z-direction for cylindrical coordinates; positive = OUT of model domain. A matrix of dimension Nr*Nteta [M/L2/T].

flux.phi.up

flux across the upstream boundary in phi-direction, for spherical coordinates; positive = OUT of model domain. A matrix of dimension Nr*Nteta [M/L2/T].

flux.phi.down

flux across the downstream boundary in phi-direction, for spherical coordinates; positive = OUT of model domain. A matrix of dimension Nr*Nteta [M/L2/T].

Details

tran.cylindrical performs (diffusive) transport in cylindrical coordinates

tran.spherical performs (diffusive) transport in spherical coordinates

The boundary conditions are either

• (1) zero gradient

• (2) fixed concentration

• (3) fixed flux

• (4) cyclic boundary

This is also the order of priority. The cyclic boundary overrules the other. If fixed concentration, fixed flux, and cyclicBnd are NULL then the boundary is zero-gradient

A cyclic boundary condition has concentration and flux at upstream and downstream boundary the same. It is useful mainly for the theta and phi direction.

** Do not expect too much of this equation: do not try to use it with many boxes **

See Also

tran.polar for a discretisation of 2-D transport equations in polar coordinates

tran.1D, tran.2D, tran.3D

Examples

# NOT RUN {
## =============================================================================
## Testing the functions
## =============================================================================
# Grid definition
r.N     <- 4   # number of cells in r-direction
theta.N <- 6   # number of cells in theta-direction
z.N     <- 3   # number of cells in z-direction

D       <- 100 # diffusion coefficient
 
r      <- seq(0,   8, len = r.N+1)       # cell size r-direction [cm]
theta  <- seq(0,2*pi, len = theta.N+1)   # theta-direction - theta: from 0, 2pi
phi    <- seq(0,2*pi, len = z.N+1)       # phi-direction (0,2pi)
z      <- seq(0,5, len = z.N+1)          # cell size z-direction [cm]
 
# Intial conditions 
C <- array(dim = c(r.N, theta.N, z.N), data = 0)

# Concentration boundary conditions
tran.cylindrical (C = C, D.r = D, D.theta = D, 
  C.r.up = 1, C.r.down = 1,
  C.theta.up = 1, C.theta.down = 1, 
  C.z.up = 1, C.z.down = 1,
  r = r, theta = theta, z = z )

tran.spherical (C = C, D.r = D, D.theta = D, 
  C.r.up = 1, C.r.down = 1, C.theta.up = 1, C.theta.down = 1, 
  C.phi.up = 1, C.phi.down = 1,
  r = r, theta = theta, phi = phi)

# Flux boundary conditions
tran.cylindrical(C = C, D.r = D, r = r, theta = theta, z = z,
  flux.r.up = 10, flux.r.down = 10,
  flux.theta.up = 10, flux.theta.down = 10,
  flux.z.up = 10, flux.z.down = 10)

tran.spherical(C = C, D.r = D, r = r, theta = theta, phi = phi,
  flux.r.up = 10, flux.r.down = 10,
  flux.theta.up = 10, flux.theta.down = 10,
  flux.phi.up = 10, flux.phi.down = 10)

# cyclic boundary conditions
tran.cylindrical(C = C, D.r = D, r = r, theta = theta, z = z,
  cyclicBnd = 1:3)
tran.spherical(C = C, D.r = D, r = r, theta = theta, phi = phi,
  cyclicBnd = 1:3)

# zero-gradient boundary conditions
tran.cylindrical(C = C, D.r = D, r = r, theta = theta, z = z)
tran.spherical(C = C, D.r = D, r = r, theta = theta, phi = phi)

## =============================================================================
## A model with diffusion and first-order consumption
## =============================================================================

N     <- 10          # number of grid cells
rr    <- 0.005       # consumption rate
D     <- 400

r       <- seq (2, 4, len = N+1)
theta   <- seq (0, 2*pi, len = N+1)
z       <- seq (0, 3, len = N+1)
phi     <- seq (0, 2*pi, len = N+1)

# The model equations
Diffcylin <- function (t, y, parms)  {
  CONC  <- array(dim = c(N, N, N), data = y)
  tran  <- tran.cylindrical(CONC, 
        D.r = D, D.theta = D, D.z = D,
        r = r, theta = theta, z = z,
        C.r.up = 0,  C.r.down = 1,
        cyclicBnd = 2)
  dCONC <- tran$dC  - rr * CONC
  return (list(dCONC))
}

Diffspher <- function (t, y, parms)  {
  CONC  <- array(dim = c(N, N, N), data = y)
  tran  <- tran.spherical (CONC, 
        D.r = D, D.theta = D, D.phi = D,
        r = r, theta = theta, phi = phi,
        C.r.up = 0,  C.r.down = 1,
        cyclicBnd = 2:3)
  dCONC <- tran$dC  - rr * CONC
  return (list(dCONC))
}

# initial condition: 0 everywhere, except in central point
y   <- array(dim = c(N, N, N), data = 0)
N2  <- ceiling(N/2)

y[N2, N2, N2] <- 100  # initial concentration in the central point...

# solve to steady-state; cyclicBnd = 2, 
outcyl <- steady.3D (y = y, func = Diffcylin, parms = NULL,
                  dim = c(N, N, N), lrw = 1e6, cyclicBnd = 2)

STDcyl <- array(dim = c(N, N, N), data = outcyl$y)
image(STDcyl[,,1])

# For spherical coordinates, cyclic Bnd = 2, 3
outspher <- steady.3D (y = y, func = Diffspher, parms = NULL, pos=TRUE,
                  dim = c(N, N, N), lrw = 1e6, cyclicBnd = 2:3)

#STDspher <- array(dim = c(N, N, N), data = outspher$y)
#image(STDspher[,,1])

# }
# NOT RUN {
  image(outspher)
# }