xenial (1) r.solute.transport.1grass.gz

Provided by: grass-doc_7.0.3-1build1_all bug

NAME

       r.solute.transport   -  Numerical  calculation  program  for  transient,  confined  and unconfined solute
       transport in two dimensions

KEYWORDS

       raster, hydrology, solute transport

SYNOPSIS

       r.solute.transport
       r.solute.transport --help
       r.solute.transport [-fc]  c=name  phead=name  hc_x=name  hc_y=name  status=name  diff_x=name  diff_y=name
       [q=name]    [cin=name]   cs=name  rd=name nf=name top=name bottom=name output=name  [vx=name]   [vy=name]
       dtime=float  [maxit=integer]   [error=float]   [solver=name]    [relax=float]    [al=float]    [at=float]
       [loops=float]   [stab=string]   [--overwrite]  [--help]  [--verbose]  [--quiet]  [--ui]

   Flags:
       -f
           Use a full filled quadratic linear equation system, default is a sparse linear equation system.

       -c
           Use the Courant-Friedrichs-Lewy criteria for time step calculation

       --overwrite
           Allow output files to overwrite existing files

       --help
           Print usage summary

       --verbose
           Verbose module output

       --quiet
           Quiet module output

       --ui
           Force launching GUI dialog

   Parameters:
       c=name [required]
           The initial concentration in [kg/m^3]

       phead=name [required]
           The piezometric head in [m]

       hc_x=name [required]
           The x-part of the hydraulic conductivity tensor in [m/s]

       hc_y=name [required]
           The y-part of the hydraulic conductivity tensor in [m/s]

       status=name [required]
           The  status  for  each  cell,  =  0 - inactive cell, 1 - active cell, 2 - dirichlet- and 3 - transfer
           boundary condition

       diff_x=name [required]
           The x-part of the diffusion tensor in [m^2/s]

       diff_y=name [required]
           The y-part of the diffusion tensor in [m^2/s]

       q=name
           Groundwater sources and sinks in [m^3/s]

       cin=name
           Concentration sources and sinks bounded to a water source or sink in [kg/s]

       cs=name [required]
           Concentration of inner sources and inner sinks in [kg/s] (i.e. a chemical reaction)

       rd=name [required]
           Retardation factor [-]

       nf=name [required]
           Effective porosity [-]

       top=name [required]
           Top surface of the aquifer in [m]

       bottom=name [required]
           Bottom surface of the aquifer in [m]

       output=name [required]
           The resulting concentration of the numerical solute transport calculation will  be  written  to  this
           map. [kg/m^3]

       vx=name
           Calculate and store the groundwater filter velocity vector part in x direction [m/s]

       vy=name
           Calculate and store the groundwater filter velocity vector part in y direction [m/s]

       dtime=float [required]
           The calculation time in seconds
           Default: 86400

       maxit=integer
           Maximum number of iteration used to solve the linear equation system
           Default: 10000

       error=float
           Error break criteria for iterative solver
           Default: 0.000001

       solver=name
           The type of solver which should solve the linear equation system
           Options: gauss, lu, jacobi, sor, bicgstab
           Default: bicgstab

       relax=float
           The relaxation parameter used by the jacobi and sor solver for speedup or stabilizing
           Default: 1

       al=float
           The longditudinal dispersivity length. [m]
           Default: 0.0

       at=float
           The transversal dispersivity length. [m]
           Default: 0.0

       loops=float
           Use this number of time loops if the CFL flag is off. The timestep will become dt/loops.
           Default: 1

       stab=string
           Set the flow stabilizing scheme (full or exponential upwinding).
           Options: full, exp
           Default: full

DESCRIPTION

       This  numerical  program  calculates  numerical  implicit  transient and steady state solute transport in
       porous media in the saturated zone of an aquifer. The computation is based on raster maps and the current
       region  settings.  All  initial- and boundary-conditions must be provided as raster maps. The unit in the
       location must be meters.

       This module is sensitive to mask settings. All cells which are outside the mask are ignored  and  handled
       as no flow boundaries.
       This  module  calculates  the concentration of the solution and optional the velocity field, based on the
       hydraulic conductivity, the effective porosity and the initial piecometric heads.  The vector  components
       can be visualized with paraview if they are exported with r.out.vtk.
       Use r.gwflow to compute the piezometric heights of the aquifer. The piezometric heights and the hydraulic
       conductivity are used to compute the flow direction and the mean velocity of the  groundwater.   This  is
       the base of the solute transport computation.
       The  solute  transport will always be calculated transient.  For stady state computation set the timestep
       to a large number (billions of seconds).
       To reduce the numerical dispersion, which is a consequence of the convection term and the  finite  volume
       discretization, you can use small time steps and choose between full and exponential upwinding.

NOTES

       The  solute  transport calculation is based on a diffusion/convection partial differential equation and a
       numerical implicit finite volume discretization. Specific for this kind of differential equation  is  the
       combination  of  a  diffusion/dispersion  term  and  a convection term.  The discretization results in an
       unsymmetric linear equation system in form of Ax = b, which must be solved. The solute transport  partial
       differential equation is of the following form:

       (dc/dt)*R = div ( D grad c - uc) + cs -q/nf(c - c_in)

           •   c -- the concentration [kg/m^3]

           •   u -- vector of mean groundwater flow velocity

           •   dt -- the time step for transient calculation in seconds [s]

           •   R -- the linear retardation coefficient [-]

           •   D -- the diffusion and dispersion tensor [m^2/s]

           •   cs -- inner concentration sources/sinks [kg/m^3]

           •   c_in -- the solute concentration of influent water [kg/m^3]

           •   q -- inner well sources/sinks [m^3/s]

           •   nf -- the effective porosity [-]
       Three  different boundary conditions are implemented, the Dirichlet, Transmission and Neumann conditions.
       The calculation and boundary status of single cells can be set with the status map.  The following states
       are supportet:

           •   0  ==  inactive  - the cell with status 0 will not be calulated, active cells will have a no flow
               boundary to an inactive cell

           •   1 == active - this cell is used for sloute transport calculation, inner sources  can  be  defined
               for those cells

           •   2  ==  Dirichlet  -  cells of this type will have a fixed concentration value which do not change
               over time

           •   3 == Transmission - cells of this type should be placed on out-flow boundaries to assure the flow
               of the solute stream out
       Note  that  all  required  raster maps are read into main memory. Additionally the linear equation system
       will be allocated, so the memory consumption of this module rapidely grow with  the  size  of  the  input
       maps.
       The  resulting  linear  equation  system  Ax  =  b can be solved with several solvers.  Several iterative
       solvers with unsymmetric sparse and quadratic matrices support are implemented.  The jacobi  method,  the
       Gauss-Seidel  method  and  the biconjugate gradients-stabilized (bicgstab) method.  Additionally a direct
       Gauss solver and LU solver are available. Those direct solvers only work with quadratic matrices,  so  be
       careful  using  them with large maps (maps of size 10.000 cells will need more than one gigabyte of ram).
       Always prefer a sparse matrix solver.

EXAMPLE

       Use this small python script to create a working groundwater flow / solute transport area and data.  Make
       sure you are not in a lat/lon projection.
       #!/usr/bin/env python
       # This is an example script how groundwater flow and solute transport are
       # computed within grass
       import sys
       import os
       import grass.script as grass
       # Overwrite existing maps
       grass.run_command("g.gisenv", set="OVERWRITE=1")
       grass.message(_("Set the region"))
       # The area is 200m x 100m with a cell size of 1m x 1m
       grass.run_command("g.region", res=1, res3=1, t=10, b=0, n=100, s=0, w=0, e=200)
       grass.run_command("r.mapcalc", expression="phead = if(col() == 1 , 50, 40)")
       grass.run_command("r.mapcalc", expression="phead = if(col() ==200  , 45 + row()/40, phead)")
       grass.run_command("r.mapcalc", expression="status = if(col() == 1 || col() == 200 , 2, 1)")
       grass.run_command("r.mapcalc", expression="well = if((row() == 50 && col() == 175) || (row() == 10 && col() == 135) , -0.001, 0)")
       grass.run_command("r.mapcalc", expression="hydcond = 0.00005")
       grass.run_command("r.mapcalc", expression="recharge = 0")
       grass.run_command("r.mapcalc", expression="top_conf = 20")
       grass.run_command("r.mapcalc", expression="bottom = 0")
       grass.run_command("r.mapcalc", expression="poros = 0.17")
       grass.run_command("r.mapcalc", expression="syield = 0.0001")
       grass.run_command("r.mapcalc", expression="null = 0.0")
       grass.message(_("Compute a steady state groundwater flow"))
       grass.run_command("r.gwflow", solver="cg", top="top_conf", bottom="bottom", phead="phead",\
         status="status", hc_x="hydcond", hc_y="hydcond", q="well", s="syield",\
         recharge="recharge", output="gwresult_conf", dt=8640000000000, type="confined")
       grass.message(_("generate the transport data"))
       grass.run_command("r.mapcalc", expression="c = if(col() == 15 && row() == 75 , 500.0, 0.0)")
       grass.run_command("r.mapcalc", expression="cs = if(col() == 15 && row() == 75 , 0.0, 0.0)")
       grass.run_command("r.mapcalc", expression="tstatus = if(col() == 1 || col() == 200 , 3, 1)")
       grass.run_command("r.mapcalc", expression="diff = 0.0000001")
       grass.run_command("r.mapcalc", expression="R = 1.0")
       # Compute the initial state
       grass.run_command("r.solute.transport", solver="bicgstab", top="top_conf",\
         bottom="bottom", phead="gwresult_conf", status="tstatus", hc_x="hydcond", hc_y="hydcond",\
         rd="R", cs="cs", q="well", nf="poros", output="stresult_conf_0", dt=3600, diff_x="diff",\
         diff_y="diff", c="c", al=0.1, at=0.01)
       # Compute the solute transport for 300 days in 10 day steps
       for dt in range(30):
           grass.run_command("r.solute.transport", solver="bicgstab", top="top_conf",\
           bottom="bottom", phead="gwresult_conf", status="tstatus", hc_x="hydcond", hc_y="hydcond",\
           rd="R", cs="cs", q="well", nf="poros", output="stresult_conf_" + str(dt + 1), dt=864000, diff_x="diff",\
           diff_y="diff", c="stresult_conf_" + str(dt), al=0.1, at=0.01)

SEE ALSO

       r.gwflow
       r3.gwflow
       r.out.vtk

AUTHOR

       Sören Gebbert

       This  work  is based on the Diploma Thesis of Sören Gebbert available here at Technical University Berlin
       in Germany.

       Last changed: $Date: 2015-01-25 18:56:33 +0100 (Sun, 25 Jan 2015) $

       Main index | Raster index | Topics index | Keywords index | Full index

       © 2003-2016 GRASS Development Team, GRASS GIS 7.0.3 Reference Manual