Cranknicolson: Cranknicolson
In geocryology/PermafrostTools: PermafrostTools

Description Usage Arguments Details Value Author(s) Examples

Implicit solution of the 1D heat conduction equation.

1 2	Cranknicolson(Tj, kj, cj, gj, zj, dt, bcl, bcu, bcutype = "NEUMANN", bcltype = "NEUMANN")

`Tj`	Temperature [K or C], array with dimensions (Z)
`kj`	Thermal conductivity [W m-1 K-1], array with dimensions (Z)
`cj`	Volumetric heat capacity [J m-3 K-1], array with dimensions (Z)
`gj`	Volumetric heat source [W], array with dimensions (Z) This is ignored for boundary nodes with a Dirichlet condition.
`zj`	Z-levels increasing upwards [m], array with dimensions (Z)
`dt`	Time step [s]
`bcl`	Lower boundary condition as either [W m-2] (Neumann) or [C or K] (Dirichlet).
`bcu`	Upper boundary condition as either [W m-2] (Neumann) or [C or K] (Dirichlet).
`bcutype`	Type of upper boundary condition used. Can be "NEUMANN" (standard, fixed heat flux) or "DIRICHLET" (fixed temperature).
`bcltype`	Type of upper boundary condition used. Can be "NEUMANN" (standard, fixed heat flux) or "DIRICHLET" (fixed temperature).

1-dimensional finite-difference (crank-nicolson) heat conduction scheme that operates on variable depth-discretization. The boundary conditions can be selected as either Dirichlet-type (prescribed temperature) or Neumann-type (prescribed energy flux).

Returns An array (Z) of temperatures [C or K] resulting from forcing the system one time step further with the given bounday conditions.

Stephan Gruber <stephan.gruber@carleton.ca> TODO: introduce heat transfer coefficient for boundary condition

Tj <- c(  1,    1,    1,    1)
kj <- c(  2,    2,    2,    2)
cj <- c(2e6,  2e6,  2e6,  2e6)
gj <- c(  0,    0,    0,    0)
zj <- c(  0,   -1,   -2,   -3)
dt <- 30000
bcl <-  1
bcu <- -1

#set up plot
plot(Tj, zj, type="l", xlim=c(-1,2))

#loop over heat conduction and add new lines
for (n in 1:100) {
  Tj <- Cranknicolson(Tj, kj, cj, gj, zj, dt, bcl, bcu, 
                      bcutype = "DIRICHLET", bcltype = "NEUMANN")
  lines(Tj, zj)
}