[GRASS-SVN] r71103 - in grass-addons/grass7/raster/r.green/r.green.gshp: . libgshp r.green.gshp.technical

svn_grass at osgeo.org svn_grass at osgeo.org
Sun May 21 23:18:32 PDT 2017


Author: zarch
Date: 2017-05-21 23:18:32 -0700 (Sun, 21 May 2017)
New Revision: 71103

Added:
   grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/ashrae.py
   grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/
   grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/Makefile
   grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.html
   grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.py
Modified:
   grass-addons/grass7/raster/r.green/r.green.gshp/Makefile
   grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/Makefile
   grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.html
Log:
r.green: Add ASHRAE method to assess the length of the BHE

Modified: grass-addons/grass7/raster/r.green/r.green.gshp/Makefile
===================================================================
--- grass-addons/grass7/raster/r.green/r.green.gshp/Makefile	2017-05-22 06:11:21 UTC (rev 71102)
+++ grass-addons/grass7/raster/r.green/r.green.gshp/Makefile	2017-05-22 06:18:32 UTC (rev 71103)
@@ -3,7 +3,8 @@
 PGM=r.green.gshp
 
 SUBDIRS = libgshp \
-	r.green.gshp.theoretical
+	r.green.gshp.theoretical \
+	r.green.gshp.technical
 
 include $(MODULE_TOPDIR)/include/Make/Dir.make
 

Modified: grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/Makefile
===================================================================
--- grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/Makefile	2017-05-22 06:11:21 UTC (rev 71102)
+++ grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/Makefile	2017-05-22 06:18:32 UTC (rev 71103)
@@ -3,7 +3,7 @@
 include $(MODULE_TOPDIR)/include/Make/Other.make
 include $(MODULE_TOPDIR)/include/Make/Python.make
 
-MODULES = gpot __init__
+MODULES = gpot ashrae __init__
 
 PGM = r.green
 LIBDIR = libgshp

Added: grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/ashrae.py
===================================================================
--- grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/ashrae.py	                        (rev 0)
+++ grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/ashrae.py	2017-05-22 06:18:32 UTC (rev 71103)
@@ -0,0 +1,1216 @@
+#!/usr/bin/env python3
+# -*- coding: utf-8 -*-
+"""
+ASHREA methodology
+------------------
+
+"""
+from collections import namedtuple
+from os import getpid
+
+from numpy import log, pi
+
+# import grass libraries
+from grass.script import raster as grast
+from grass.script import core as gcore
+
+from grass.pygrass.utils import get_mapset_raster
+
+
+BASENAME = 'tmprgreen{pid:05d}_'.format(pid=getpid())
+
+# Coefficients for Tp correlation
+TP = [7.8189e+00, -6.4270e+01, 1.5387e+02, -8.4809e+01, 3.4610e+00,
+      -9.4753e-01, -6.0416e-02, 1.5631e+00, -8.9416e-03, 1.9061e-05,
+      -2.2890e+00, 1.0187e-01, 6.5690e-03, -4.0918e+01, 1.5557e+01,
+      -1.9107e+01, 1.0529e-01, 2.5501e+01, -2.1177e+00, 7.7529e+01,
+      -5.0454e+01, 7.6352e+01, -5.3719e-01, -1.3200e+02, 1.2878e+01,
+      1.2697e-01, -4.0284e-04, -7.2065e-02, 9.5184e-04, -2.4167e-02,
+      9.6811e-05, 2.8317e-02, -1.0905e-03, 1.2207e-01, -7.1050e-03,
+      -1.1129e-03, -4.5566e-04]
+
+
+FUNCTIONS = {'6h': [0.6619352, -4.815693, 15.03571, -0.09879421, 0.02917889,
+                    0.1138498, 0.005610933, 0.7796329, -0.3243880,
+                    -0.01824101],
+             '1m': [0.4132728, 0.2912981, 0.07589286, 0.1563978, -0.2289355,
+                    -0.004927554, -0.002694979, -0.6380360, 0.2950815,
+                    0.1493320],
+             '10y': [0.3057646, 0.08987446, -0.09151786, -0.03872451,
+                     0.1690853, -0.02881681, -0.002886584, -0.1723169,
+                     0.03112034, -0.1188438]}
+
+
+# define named tuple to struct data/characteristics
+GroundLoads = namedtuple('GroundLoads',
+                         field_names=['hourly', 'monthly', 'yearly'])
+
+
+GroundProperties = namedtuple('GroundProperties',
+                              field_names=['conductivity', 'diffusivity',
+
+                                           'temperature'])
+FluidProperties = namedtuple('FluidProperties',
+                             field_names=['capacity', 'massflow', 'inlettemp'])
+
+
+Borehole = namedtuple('Borehole',
+                      field_names=['radius',
+                                   'pipe_inner_radius', 'pipe_outer_radius',
+                                   'k_grout', 'k_pipe',
+                                   'distance', 'convection'])
+
+
+BoreholeExchanger = namedtuple('BoreholeExchanger',
+                               field_names=['ground_loads',
+                                            'ground',
+                                            'fluid',
+                                            'borehole'])
+
+BoreholeField = namedtuple('BoreholeField',
+                           field_names=['distance', 'number', 'ratio', 'bhe'])
+
+
+InfoVars = namedtuple('InfoVars', field_names=['long_term', 'medium_term',
+                                               'short_term', 'fluid_temp',
+                                               'resistence'])
+
+
+def exists(raster_name):
+    """Return True if the map already exist
+    """
+    raster, mapset = (raster_name.split('@') if '@' in raster_name else
+                      (raster_name, ''))
+    return bool(get_mapset_raster(raster_name, mapset))
+
+
+def rename(mtype, from_mname, to_mname):
+    from_to = '{from_name},{to_name}'
+    gcore.run_command('g.rename',
+                      **{mtype: from_to.format(from_name=from_mname,
+                                               to_name=to_mname)})
+
+
+def ground_resistence(ground, borehole, period='6h'):
+    """Return the effective ground thermal resistances:
+
+    :math:`R_{period}` [m K W-1]
+
+    Parameters
+    ----------
+
+    ground: GroundProperties
+        GroundProperties tuple with the main characteristics of the ground
+    borehole: rbore [m]
+        borehole radius
+    period: str
+        Six periods are supported:
+            * 6h (short term),
+            * 1m (medium term),
+            * 10y (long term)
+
+    Example
+    -------
+    >>> ground = GroundProperties(conductivity=2, diffusivity=0.086,
+    ...                           temperature=10.)
+    >>> ground_resistence(ground, borehole=0.06, period='6h')
+    0.11445903013141503
+    >>> ground_resistence(ground, borehole=0.06, period='1m')
+    0.18020547706451287
+    >>> ground_resistence(ground, borehole=0.06, period='10y')
+    0.19083864808240175
+
+    >>> ground = GroundProperties(conductivity=2.25, diffusivity=0.06757039008,
+    ...                           temperature=10.)
+    >>> ground_resistence(ground, borehole=0.054, period='6h')
+    0.10067477057511355
+    >>> ground_resistence(ground, borehole=0.054, period='1m')
+    0.1600011827001955
+    >>> ground_resistence(ground, borehole=0.054, period='10y')
+    0.16963475237210959
+    """
+    c0, c1, c2, c3, c4, c5, c6, c7, c8, c9 = FUNCTIONS[period]
+    conductivity, diffusivity, temperature = ground
+    log_diff = log(diffusivity)
+    return (1. / conductivity *
+            (c0 +
+             c1 * borehole +
+             c2 * borehole**2. +
+             c3 * diffusivity +
+             c4 * diffusivity**2. +
+             c5 * log_diff +
+             c6 * log_diff**2 +
+             c7 * borehole * diffusivity +
+             c8 * borehole * log_diff +
+             c9 * diffusivity * log_diff))
+
+
+def _log(out, value, execute=True, show=False, **kwargs):
+    """Return the natural log of a number/raster
+
+    Example
+    -------
+    >>> _log('log_elev', 'elevation', execute=False, show=True)
+    log_elev = log(elevation)
+    'log_elev'
+    >>> _log('log_elev', 1, execute=False)
+    0.0
+    """
+    if isinstance(value, str):
+        # use rmapcalc and return raster name
+        rcmd = "{out} = log({name})".format(out=out, name=value)
+        if show:
+            print(rcmd)
+        if execute:
+            grast.mapcalc(rcmd, **kwargs)
+        return out
+    else:
+        return log(value)
+
+
+def r_ground_resistence(out, ground, borehole, period='6h',
+                        execute=True, **kwargs):
+    """Return the effective ground thermal resistances:
+
+    :math:`R_{period}` [m K W-1]
+
+    Parameters
+    ----------
+
+    ground: GroundProperties
+        GroundProperties tuple with the main characteristics of the ground
+    borehole: rbore [m]
+        borehole radius
+    period: str
+        Three periods are supported:
+            * 6h (short term),
+            * 1m (medium term),
+            * 10y (long term)
+
+    Example
+    -------
+    >>> ground = GroundProperties(conductivity=2, diffusivity=0.086,
+    ...                           temperature=10.)
+    >>> r_ground_resistence('g_resistence_6h', ground, borehole=0.06,
+    ...                     period='6h', execute=False)
+    ...                                       # doctest: +NORMALIZE_WHITESPACE
+    'g_resistence_6h = (1. / 2 *
+                        (0.6619352 +
+                         -4.815693 * 0.06 +
+                         15.03571 * 0.06^2. +
+                         -0.09879421 * 0.086 +
+                         0.02917889 * 0.086^2. +
+                         0.1138498 * -2.45340798273 +
+                         0.005610933 * -2.45340798273^2 +
+                         0.7796329 * 0.06 * 0.086 +
+                         -0.324388 * 0.06 * -2.45340798273 +
+                         -0.01824101 * 0.086 * -2.45340798273))'
+    >>> r_ground_resistence('g_resistence_1m', ground, borehole=0.06,
+    ...                     period='1m', execute=False)
+    ...                                       # doctest: +NORMALIZE_WHITESPACE
+    'g_resistence_1m = (1. / 2 *
+                        (0.4132728 +
+                         0.2912981 * 0.06 +
+                         0.07589286 * 0.06^2. +
+                         0.1563978 * 0.086 +
+                         -0.2289355 * 0.086^2. +
+                         -0.004927554 * -2.45340798273 +
+                         -0.002694979 * -2.45340798273^2 +
+                         -0.638036 * 0.06 * 0.086 +
+                         0.2950815 * 0.06 * -2.45340798273 +
+                         0.149332 * 0.086 * -2.45340798273))'
+    >>> r_ground_resistence('g_resistence_10y', ground, borehole=0.06,
+    ...                     period='10y', execute=False)
+    ...                                       # doctest: +NORMALIZE_WHITESPACE
+    'g_resistence_10y = (1. / 2 *
+                         (0.3057646 +
+                          0.08987446 * 0.06 +
+                          -0.09151786 * 0.06^2. +
+                          -0.03872451 * 0.086 +
+                          0.1690853 * 0.086^2. +
+                          -0.02881681 * -2.45340798273 +
+                          -0.002886584 * -2.45340798273^2 +
+                          -0.1723169 * 0.06 * 0.086 +
+                          0.03112034 * 0.06 * -2.45340798273 +
+                          -0.1188438 * 0.086 * -2.45340798273))'
+
+
+    >>> ground = GroundProperties(conductivity='g_conductivity',
+    ...                           diffusivity='g_diffusivity',
+    ...                           temperature='g_temperature')
+    >>> r_ground_resistence('g_resistence_6h', ground, borehole=0.06,
+    ...                     period='6h', execute=False)
+    ...                                       # doctest: +NORMALIZE_WHITESPACE
+    'g_resistence_6h = (1. / g_conductivity *
+                        (0.6619352 +
+                         -4.815693 * 0.06 +
+                         15.03571 * 0.06^2. +
+                         -0.09879421 * g_diffusivity +
+                         0.02917889 * g_diffusivity^2. +
+                         0.1138498 * g_resistence_6h_log +
+                         0.005610933 * g_resistence_6h_log^2 +
+                         0.7796329 * 0.06 * g_diffusivity +
+                         -0.324388 * 0.06 * g_resistence_6h_log +
+                         -0.01824101 * g_diffusivity * g_resistence_6h_log))'
+    >>> r_ground_resistence('g_resistence_1m', ground, borehole=0.06,
+    ...                     period='1m', execute=False)
+    ...                                       # doctest: +NORMALIZE_WHITESPACE
+    'g_resistence_1m = (1. / g_conductivity *
+                        (0.4132728 +
+                         0.2912981 * 0.06 +
+                         0.07589286 * 0.06^2. +
+                         0.1563978 * g_diffusivity +
+                         -0.2289355 * g_diffusivity^2. +
+                         -0.004927554 * g_resistence_1m_log +
+                         -0.002694979 * g_resistence_1m_log^2 +
+                         -0.638036 * 0.06 * g_diffusivity +
+                         0.2950815 * 0.06 * g_resistence_1m_log +
+                         0.149332 * g_diffusivity * g_resistence_1m_log))'
+    >>> r_ground_resistence('g_resistence_10y', ground, borehole=0.06,
+    ...                     period='10y', execute=False)
+    ...                                       # doctest: +NORMALIZE_WHITESPACE
+    'g_resistence_10y = (1. / g_conductivity *
+                         (0.3057646 +
+                          0.08987446 * 0.06 +
+                          -0.09151786 * 0.06^2. +
+                          -0.03872451 * g_diffusivity +
+                          0.1690853 * g_diffusivity^2. +
+                          -0.02881681 * g_resistence_10y_log +
+                          -0.002886584 * g_resistence_10y_log^2 +
+                          -0.1723169 * 0.06 * g_diffusivity +
+                          0.03112034 * 0.06 * g_resistence_10y_log +
+                          -0.1188438 * g_diffusivity * g_resistence_10y_log))'
+    """
+    c0, c1, c2, c3, c4, c5, c6, c7, c8, c9 = FUNCTIONS[period]
+    conductivity, diffusivity, temperature = ground
+    log_diff = _log(out + '_log', diffusivity, execute=execute, **kwargs)
+    res = ("{out} = (1. / {conductivity} *"
+           " ({c0} +"
+           "  {c1} * {borehole} +"
+           "  {c2} * {borehole}^2. +"
+           "  {c3} * {diffusivity} +"
+           "  {c4} * {diffusivity}^2. +"
+           "  {c5} * {log_diff} +"
+           "  {c6} * {log_diff}^2 +"
+           "  {c7} * {borehole} * {diffusivity} +"
+           "  {c8} * {borehole} * {log_diff} +"
+           "  {c9} * {diffusivity} * {log_diff}))")
+    rcmd = res.format(out=out, c0=c0, c1=c1, c2=c2, c3=c3, c4=c4, c5=c5,
+                      c6=c6, c7=c7, c8=c8, c9=c9, conductivity=conductivity,
+                      diffusivity=diffusivity, borehole=borehole,
+                      log_diff=log_diff)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def fluid_temperature_hp_outlet(fluid, peak):
+    """Return heat pump outlet temperature :math:`T_{outHP}` [°C]
+
+    Parameters
+    ----------
+    fluid: FluidProperties
+        FluidProperties named tuple with the  fluid characteristics
+    peak: qh [W]
+        peak hourly ground load
+
+    Example
+    -------
+    >>> fluid = FluidProperties(capacity=4200, massflow=0.050,
+    ...                         inlettemp=40.2)
+    >>> fluid_temperature_hp_outlet(fluid, 12000)
+    44.96190476190476
+
+    >>> fluid = FluidProperties(capacity=4000, massflow=0.074,
+    ...                         inlettemp=4.44)
+    >>> fluid_temperature_hp_outlet(fluid, -392250)
+    1.0616216216216219
+    """
+    return (fluid.inlettemp +
+            peak / (fluid.massflow * abs(peak) / 1000. * fluid.capacity))
+
+
+def r_fluid_temperature_hp_outlet(out, fluid, peak, execute=True, **kwargs):
+    """Return heat pump outlet temperature :math:`T_{outHP}` [°C]
+
+    Parameters
+    ----------
+    fluid: FluidProperties
+        FluidProperties named tuple with the  fluid characteristics
+    peak: qh [W]
+        peak hourly ground load
+
+    Example
+    -------
+    >>> fluid = FluidProperties(capacity=4200, massflow=0.050,
+    ...                         inlettemp=40.2)
+    >>> r_fluid_temperature_hp_outlet('fluid_tmp_hp_outlet', fluid, 'peak',
+    ...                               execute=False)
+    'fluid_tmp_hp_outlet = 40.2 + peak / (0.05 * abs(peak) / 1000. * 4200)'
+    """
+    res = ("{out} = {fluid_inlettemp} + "
+           "{peak} / "
+           "({fluid_massflow} * abs({peak}) / 1000. * {fluid_capacity})")
+    rcmd = res.format(out=out, peak=peak,
+                      fluid_inlettemp=fluid.inlettemp,
+                      fluid_massflow=fluid.massflow,
+                      fluid_capacity=fluid.capacity)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def fluid_temperature_borehole(fluid, peak):
+    """Return the average fluid temperature in the borehole:
+    :math:`T_m` [°C]
+
+    Parameters
+    ----------
+    peak: qh [W]
+        peak hourly ground load
+
+    Example
+    -------
+    >>> fluid = FluidProperties(capacity=4200, massflow=0.050,
+    ...                         inlettemp=40.2)
+    >>> fluid_temperature_borehole(fluid, 12000)
+    42.58095238095238
+
+    >>> fluid = FluidProperties(capacity=4000, massflow=0.074,
+    ...                         inlettemp=4.44)
+    >>> fluid_temperature_borehole(fluid, -392250)
+    2.750810810810811
+    """
+    return (fluid.inlettemp + fluid_temperature_hp_outlet(fluid, peak)) / 2.
+
+
+def r_fluid_temperature_borehole(out, fluid, peak, execute=True, **kwargs):
+    """Return the average fluid temperature in the borehole:
+    :math:`T_m` [°C]
+
+    Parameters
+    ----------
+    peak: qh [W]
+        peak hourly ground load
+
+    Example
+    -------
+    >>> fluid = FluidProperties(capacity=4200, massflow=0.050,
+    ...                         inlettemp=40.2)
+    >>> r_fluid_temperature_borehole('fluid_temp', fluid, 'peak',
+    ...                              execute=False)
+    'fluid_temp = (40.2 + fluid_temp_hp_outlet) / 2.'
+    """
+    res = "{out} = ({fluid_inlettemp} + {fluid_temperature_hp_outlet}) / 2."
+    fluid_temperature_hp_outlet = out + '_hp_outlet'
+    r_fluid_temperature_hp_outlet(fluid_temperature_hp_outlet, fluid, peak,
+                                  execute=execute)
+    rcmd = res.format(out=out, fluid_inlettemp=fluid.inlettemp,
+                      fluid_temperature_hp_outlet=fluid_temperature_hp_outlet)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def bh_resistence_convetive(bh):
+    """Return the conductie resistence: :math:`R_{conv}` [m K W-1]
+
+    .. math::
+
+        R_{conv} = \frac{1}{2\\pi \\dot r_{pin} \\dot h_{conv}}
+
+    Example
+    -------
+    >>> bh = Borehole(radius=0.06,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...               convection=1000.)
+    >>> bh_resistence_convetive(bh)
+    0.011659702790615043
+    """
+    return 1. / (2. * pi * bh.pipe_inner_radius * bh.convection)
+
+
+def bh_resistence_pipe(bh):
+    """Return pipe resistence: :math:`R_{pipe}` [m K W-1]
+
+    .. math::
+
+        R_{pipe} = \frac{\frac{r_{pext}}{r_{pint}}}{2\\pi \\dot k_{pipe}}
+
+    Example
+    -------
+    >>> bh = Borehole(radius=0.06,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...               convection=1000.)
+    >>> bh_resistence_pipe(bh)
+    0.076420594518887289
+
+    >>> bh = Borehole(radius=0.054,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...               convection=1000.)
+    >>> bh_resistence_pipe(bh)
+    0.071325888217628128
+    """
+    return (log(bh.pipe_outer_radius / bh.pipe_inner_radius) /
+            (2. * pi * bh.k_pipe))
+
+
+def bh_resistence_grout(bh, ground_conductivity):
+    """Return grout resistence: :math:`R_{grout}` [m K W-1]
+
+    .. math::
+
+        R_{grout} = \frac{1}{4\\pi \\dot k_{grout}} \\dot
+                    \\ln(\frac{r_{bore}}{r_{pint}}) +
+                    \\ln(\frac{r_{bore}}{LU}) +
+                    \frac{k_{grout} - k_{ground}}{k_{grout} + k_{ground}} *
+                    \\ln(\frac{r_{bore}^4}{r_{bore}^4 - (\frac{LU}{2})^4})
+
+    Example
+    -------
+    >>> bh = Borehole(radius=0.06,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...               convection=1000.)
+    >>> bh_resistence_grout(bh, ground_conductivity=2.)
+    0.076114233912862317
+
+    >>> bh = Borehole(radius=0.054,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...               convection=1000.)
+    >>> bh_resistence_grout(bh, ground_conductivity=2.25)
+    0.060049831980162838
+    """
+    return (1. / (4. * pi * bh.k_grout) *
+            (log(bh.radius / bh.pipe_outer_radius) +
+             log(bh.radius / bh.distance) +
+             (bh.k_grout - ground_conductivity) /
+             (bh.k_grout + ground_conductivity) *
+             log(bh.radius**4. /
+                 (bh.radius**4. - (bh.distance / 2.)**4.))))
+
+
+def r_bh_resistence_grout(out, bh, ground_conductivity,
+                          execute=True, **kwargs):
+    """Return grout resistence: :math:`R_{grout}` [m K W-1]
+
+    .. math::
+
+        R_{grout} = \frac{1}{4\\pi \\dot k_{grout}} \\dot
+                    \\ln(\frac{r_{bore}}{r_{pint}}) +
+                    \\ln(\frac{r_{bore}}{LU}) +
+                    \frac{k_{grout} - k_{ground}}{k_{grout} + k_{ground}} *
+                    \\ln(\frac{r_{bore}^4}{r_{bore}^4 - (\frac{LU}{2})^4})
+
+    Example
+    -------
+    >>> bh = Borehole(radius=0.06,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...               convection=1000.)
+    >>> r_bh_resistence_grout('bh_resistence_grout',
+    ...                       bh, ground_conductivity='g_conductivity',
+    ...                       execute=False)   # doctest: +NORMALIZE_WHITESPACE
+    'bh_resistence_grout = (1. /
+                            (4. * 3.14159265359 * 1.5) *
+                            (log(0.06 / 0.0167) + log(0.06 / 0.0511) +
+                            (1.5 - g_conductivity) /
+                            (1.5 + g_conductivity) *
+                            log(0.06^4. /      (0.06^4. - (0.0511/ 2.)^4.))))'
+    """
+    res = ("{out} = (1. / (4. * {pi} * {bh.k_grout}) * "
+           "(log({bh.radius} / {bh.pipe_outer_radius}) +"
+           " log({bh.radius} / {bh.distance}) +"
+           " ({bh.k_grout} - {ground_conductivity}) /"
+           "  ({bh.k_grout} + {ground_conductivity}) *"
+           "  log({bh.radius}^4. /"
+           "      ({bh.radius}^4. - ({bh.distance}/ 2.)^4.))))")
+    rcmd = res.format(out=out, bh=bh, ground_conductivity=ground_conductivity,
+                      pi=pi)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def bh_resistence(bh, ground_conductivity):
+    """Return the effective borehole thermal resistance: Rb [m K W-1]
+
+    .. math::
+
+        R_{b} = R_{grout} + \frac{R_{conv} + R_{pipe}}{2}
+
+    Example
+    -------
+    >>> bh = Borehole(radius=0.06,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...               convection=1000.)
+    >>> bh_resistence(bh, ground_conductivity=2.)
+    0.12015438256761349
+
+    >>> bh = Borehole(radius=0.054,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...               convection=1000.)
+    >>> bh_resistence(bh, ground_conductivity=2.25)
+    0.10154262748428441
+    """
+    return (bh_resistence_grout(bh, ground_conductivity) +
+            (bh_resistence_convetive(bh) + bh_resistence_pipe(bh)) / 2.)
+
+
+def r_bh_resistence(out, bh, ground_conductivity, execute=True, **kwargs):
+    """Return the effective borehole thermal resistance: Rb [m K W-1]
+
+    .. math::
+
+        R_{b} = R_{grout} + \frac{R_{conv} + R_{pipe}}{2}
+
+    Example
+    -------
+    >>> bh = Borehole(radius=0.06,
+    ...               pipe_inner_radius=0.01365, pipe_outer_radius=0.0167,
+    ...               k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...               convection=1000.)
+    >>> r_bh_resistence('bh_resistence', bh,
+    ...                 ground_conductivity='g_conductivity', execute=False)
+    ...                                        # doctest: +NORMALIZE_WHITESPACE
+    'bh_resistence = (bh_resistence_grout +
+                      (0.0116597027906 + 0.0764205945189) / 2.)'
+    """
+    bh_resistence_grout = out + '_grout'
+    r_bh_resistence_grout(bh_resistence_grout, bh, ground_conductivity,
+                          execute=execute, **kwargs)
+    res = ("{out} = ({bh_resistence_grout} +"
+           "         ({bh_resistence_convetive} + {bh_resistence_pipe}) / 2.)")
+    rcmd = res.format(out=out, bh_resistence_grout=bh_resistence_grout,
+                      bh_resistence_convetive=bh_resistence_convetive(bh),
+                      bh_resistence_pipe=bh_resistence_pipe(bh))
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def bhe_length(bhe):
+    """Return the total length calculation assuming no borehole thermal
+    interference: L [m]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly=12000,
+    ...                                    monthly=6000,
+    ...                                    yearly=1500),
+    ...           ground=GroundProperties(conductivity=2, diffusivity=0.086,
+    ...                                   temperature=15.),
+    ...           fluid=FluidProperties(capacity=4200, massflow=0.050,
+    ...                                 inlettemp=40.2),
+    ...           borehole=Borehole(radius=0.06,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...                             convection=1000.)
+    ...           )
+    >>> bhe_length(bhe)
+    151.65726437306537
+
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly=-392250,
+    ...                                    monthly=-100000,
+    ...                                    yearly=-1762),
+    ...           ground=GroundProperties(conductivity=2.25, diffusivity=0.06757039008,
+    ...                                   temperature=12.41),
+    ...           fluid=FluidProperties(capacity=4000, massflow=0.074,
+    ...                                 inlettemp=4.44),
+    ...           borehole=Borehole(radius=0.054,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...                             convection=1000.)
+    ...           )
+    >>> bhe_length(bhe)
+    9899.2562646476617
+    """
+    long_term = ground_resistence(bhe.ground, bhe.borehole.radius,
+                                  period='10y')
+    medium_term = ground_resistence(bhe.ground, bhe.borehole.radius,
+                                    period='1m')
+    short_term = ground_resistence(bhe.ground, bhe.borehole.radius,
+                                   period='6h')
+    fluid_temp = fluid_temperature_borehole(bhe.fluid,
+                                            bhe.ground_loads.hourly)
+    resistence = bh_resistence(bhe.borehole, bhe.ground.conductivity)
+    return ((bhe.ground_loads.yearly * long_term +
+             bhe.ground_loads.monthly * medium_term +
+             bhe.ground_loads.hourly * short_term +
+             bhe.ground_loads.hourly * resistence) /
+            (fluid_temp - bhe.ground.temperature))
+
+
+def get_vars(out, bhe, basename, execute=True, **kwargs):
+    basename = basename if basename else BASENAME
+    # compute the temporary maps
+    long_term = basename + out + '_long_term'
+    r_ground_resistence(long_term, bhe.ground, bhe.borehole.radius,
+                        period='10y', execute=execute, **kwargs)
+    medium_term = basename + out + '_medium_term'
+    r_ground_resistence(medium_term, bhe.ground, bhe.borehole.radius,
+                        period='1m', execute=execute, **kwargs)
+    short_term = basename + out + '_short_term'
+    r_ground_resistence(short_term, bhe.ground, bhe.borehole.radius,
+                        period='6h', execute=execute, **kwargs)
+    fluid_temp = basename + out + '_fluid_temp'
+    r_fluid_temperature_borehole(fluid_temp, bhe.fluid,
+                                 bhe.ground_loads.hourly,
+                                 execute=execute, **kwargs)
+    resistence = basename + out + '_resistence'
+    r_bh_resistence(resistence, bhe.borehole, bhe.ground.conductivity,
+                    execute=execute, **kwargs)
+    return InfoVars(long_term, medium_term, short_term, fluid_temp, resistence)
+
+
+def r_bhe_length(out, bhe, infovars, execute=True, **kwargs):
+    """Return the total length calculation assuming no borehole thermal
+    interference: L [m]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly='g_loads_6h',
+    ...                                    monthly='g_loads_1m',
+    ...                                    yearly='g_loads_1y'),
+    ...           ground=GroundProperties(conductivity='g_conductivity',
+    ...                                   diffusivity='g_diffusivity',
+    ...                                   temperature='g_temperature'),
+    ...           fluid=FluidProperties(capacity=4200, massflow=0.050,
+    ...                                 inlettemp=40.2),
+    ...           borehole=Borehole(radius=0.06,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...                             convection=1000.)
+    ...           )
+    >>> infovars = InfoVars('l_term', 'm_term', 's_term', 'f_temp', 'res')
+    >>> r_bhe_length('bhe_length', bhe, infovars, execute=False)
+    ...                                        # doctest: +NORMALIZE_WHITESPACE
+    'bhe_length = ((g_loads_1y * l_term +
+                    g_loads_1m * m_term +
+                    g_loads_6h * s_term +
+                    g_loads_6h * res) /
+                   (f_temp - g_temperature))'
+    """
+    # compute the BHE length
+    res = ("{out} = (({bhe.ground_loads.yearly} * {iv.long_term} +"
+           "          {bhe.ground_loads.monthly} * {iv.medium_term} +"
+           "          {bhe.ground_loads.hourly} * {iv.short_term} +"
+           "          {bhe.ground_loads.hourly} * {iv.resistence}) /"
+           "           ({iv.fluid_temp} - {bhe.ground.temperature}))")
+    rcmd = res.format(out=out, bhe=bhe, iv=infovars)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def distance_depth_ratio(field, length):
+    """Return the distance-depth ratio: B/H
+
+    Example
+    -------
+
+    >>> field = BoreholeField(distance=6.1, number=120, ratio=1.2, bhe=None)
+    >>> distance_depth_ratio(field, length=9899.2562646476617)
+    ...                                                    # doctest: +ELLIPSIS
+    0.07394494903764...
+    """
+    return field.distance / (length / field.number)
+
+
+def r_distance_depth_ratio(out, field, length, execute=True, **kwargs):
+    """Return the distance-depth ratio: B/H
+
+    Example
+    -------
+
+    >>> field = BoreholeField(distance='f_dist', number='f_numb',
+    ...                       ratio='f_ratio', bhe=None)
+    >>> r_distance_depth_ratio('dist_depth_ratio', field,
+    ...                        length='bh_length', execute=False)
+    'dist_depth_ratio = f_dist / (bh_length / f_numb)'
+    """
+    res = "{out} = {field.distance} / ({length} / {field.number})"
+    rcmd = res.format(out=out, field=field, length=length)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def log_dimless_time(field, length):
+    """Return the logarithm of dimensionless time: [ln(t10y/ts)]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=None,
+    ...           ground=GroundProperties(conductivity=None,
+    ...                                   diffusivity=0.06757039008,
+    ...                                   temperature=None),
+    ...           fluid=None,
+    ...           borehole=None)
+    >>> field = BoreholeField(distance=6.1, number=120, ratio=1.2, bhe=bhe)
+    >>> log_dimless_time(field, length=9899.2562646476617)
+    -1.1196400189685967
+    """
+    return log(365.25 * 10 /
+               ((length / field.number)**2 /
+                (9 * field.bhe.ground.diffusivity)))
+
+
+def r_log_dimless_time(out, field, length, execute=True, **kwargs):
+    """Return the logarithm of dimensionless time: [ln(t10y/ts)]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=None,
+    ...           ground=GroundProperties(conductivity=None,
+    ...                                   diffusivity='g_diffusivity',
+    ...                                   temperature=None),
+    ...           fluid=None,
+    ...           borehole=None)
+    >>> field = BoreholeField(distance='f_dist', number='f_numb',
+    ...                       ratio='f_ratio', bhe=bhe)
+    >>> r_log_dimless_time('ln_dimless_time', field, length='bh_length',
+    ...                    execute=False)      # doctest: +NORMALIZE_WHITESPACE
+    'ln_dimless_time = log(365.25 * 10 /
+                           ((bh_length / f_numb)^2 / (9 * g_diffusivity)))'
+    """
+    res = ("{out} = log(365.25 * 10 /"
+           "            (({length} / {field.number})^2 /"
+           "             (9 * {field.bhe.ground.diffusivity})))")
+    rcmd = res.format(out=out, field=field, length=length)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def _temperature_penality(yearly_load, ground_conductivity, length,
+                          dd_ratio, ln_dim_time, number, ratio):
+    """Return the temperature penality of a BHE field
+
+    Example
+    -------
+    >>> _temperature_penality(yearly_load=-1762,
+    ...                       ground_conductivity=2.25,
+    ...                       length=9899.25626464766,
+    ...                       dd_ratio=0.0739449490376491,
+    ...                       ln_dim_time=-1.1196400189686,
+    ...                       number=120, ratio=1.2)
+    -0.24002529024227007
+    """
+    return (yearly_load /
+            (2 * pi * ground_conductivity * length) *
+            (TP[0] +
+             TP[1] * dd_ratio +
+             TP[2] * dd_ratio**2 +
+             TP[3] * dd_ratio**3 +
+             TP[4] * ln_dim_time +
+             TP[5] * ln_dim_time**2 +
+             TP[6] * ln_dim_time**3 +
+             TP[7] * number +
+             TP[8] * number**2 +
+             TP[9] * number**3 +
+             TP[10] * ratio +
+             TP[11] * ratio**2 +
+             TP[12] * ratio**3 +
+             TP[13] * dd_ratio * ln_dim_time +
+             TP[14] * dd_ratio * ln_dim_time**2 +
+             TP[15] * dd_ratio * number +
+             TP[16] * dd_ratio * number**2 +
+             TP[17] * dd_ratio * ratio +
+             TP[18] * dd_ratio * ratio**2 +
+             TP[19] * dd_ratio**2 * ln_dim_time +
+             TP[20] * dd_ratio**2 * ln_dim_time**2 +
+             TP[21] * dd_ratio**2 * number +
+             TP[22] * dd_ratio**2 * number**2 +
+             TP[23] * dd_ratio**2 * ratio +
+             TP[24] * dd_ratio**2 * ratio**2 +
+             TP[25] * ln_dim_time * number +
+             TP[26] * ln_dim_time * number**2 +
+             TP[27] * ln_dim_time * ratio +
+             TP[28] * ln_dim_time * ratio**2 +
+             TP[29] * ln_dim_time**2 * number +
+             TP[30] * ln_dim_time**2 * number**2 +
+             TP[31] * ln_dim_time**2 * ratio +
+             TP[32] * ln_dim_time**2 * ratio**2 +
+             TP[33] * number * ratio +
+             TP[34] * number * ratio**2 +
+             TP[35] * number**2 * ratio +
+             TP[36] * number**2 * ratio**2))
+
+
+def _r_temperature_penality(out, yearly_load, ground_conductivity, length,
+                            dd_ratio, ln_dim_time, number, ratio,
+                            execute=True, **kwargs):
+    """Return the temperature penality of a BHE field
+
+    Example
+    -------
+    >>> _r_temperature_penality('temp_pen', yearly_load='y_load',
+    ...                         ground_conductivity='g_cond', length='length',
+    ...                         dd_ratio='dd_ratio', ln_dim_time='ln_dim_time',
+    ...                         number='numb', ratio='ratio', execute=False)
+    ...                                        # doctest: +NORMALIZE_WHITESPACE
+    'temp_pen = (y_load /
+                 (2 * 3.14159265359 * g_cond * length) *
+                 (7.8189 +
+                  -64.27 * dd_ratio +
+                  153.87 * dd_ratio^2 +
+                  -84.809 * dd_ratio^3 +
+                  3.461 * ln_dim_time +
+                  -0.94753 * ln_dim_time^2 +
+                  -0.060416 * ln_dim_time^3 +
+                  1.5631 * numb +
+                  -0.0089416 * numb^2 +
+                  1.9061e-05 * numb^3 +
+                  -2.289 * ratio +
+                  0.10187 * ratio^2 +
+                  0.006569 * ratio^3 +
+                  -40.918 * dd_ratio * ln_dim_time +
+                  15.557 * dd_ratio * ln_dim_time^2 +
+                  -19.107 * dd_ratio * numb +
+                  0.10529 * dd_ratio * numb^2 +
+                  25.501 * dd_ratio * ratio +
+                  -2.1177 * dd_ratio * ratio^2 +
+                  77.529 * dd_ratio^2 * ln_dim_time +
+                  -50.454 * dd_ratio^2 * ln_dim_time^2 +
+                  76.352 * dd_ratio^2 * numb +
+                  -0.53719 * dd_ratio^2 * numb^2 +
+                  -132.0 * dd_ratio^2 * ratio +
+                  12.878 * dd_ratio^2 * ratio^2 +
+                  0.12697 * ln_dim_time * numb +
+                  -0.00040284 * ln_dim_time * numb^2 +
+                  -0.072065 * ln_dim_time * ratio +
+                  0.00095184 * ln_dim_time * ratio^2 +
+                  -0.024167 * ln_dim_time^2 * numb +
+                  9.6811e-05 * ln_dim_time^2 * numb^2 +
+                  0.028317 * ln_dim_time^2 * ratio +
+                  -0.0010905 * ln_dim_time^2 * ratio^2 +
+                  0.12207 * numb * ratio +
+                  -0.007105 * numb * ratio^2 +
+                  -0.0011129 * numb^2 * ratio +
+                  -0.00045566 * numb^2 * ratio^2))'
+    """
+    res = ("{out} = ({yearly_load} /"
+           " (2 * {pi} * {ground_conductivity} * {length}) *"
+           " ({TP[0]} +"
+           "  {TP[1]} * {dd_ratio} +"
+           "  {TP[2]} * {dd_ratio}^2 +"
+           "  {TP[3]} * {dd_ratio}^3 +"
+           "  {TP[4]} * {ln_dim_time} +"
+           "  {TP[5]} * {ln_dim_time}^2 +"
+           "  {TP[6]} * {ln_dim_time}^3 +"
+           "  {TP[7]} * {number} +"
+           "  {TP[8]} * {number}^2 +"
+           "  {TP[9]} * {number}^3 +"
+           "  {TP[10]} * {ratio} +"
+           "  {TP[11]} * {ratio}^2 +"
+           "  {TP[12]} * {ratio}^3 +"
+           "  {TP[13]} * {dd_ratio} * {ln_dim_time} +"
+           "  {TP[14]} * {dd_ratio} * {ln_dim_time}^2 +"
+           "  {TP[15]} * {dd_ratio} * {number} +"
+           "  {TP[16]} * {dd_ratio} * {number}^2 +"
+           "  {TP[17]} * {dd_ratio} * {ratio} +"
+           "  {TP[18]} * {dd_ratio} * {ratio}^2 +"
+           "  {TP[19]} * {dd_ratio}^2 * {ln_dim_time} +"
+           "  {TP[20]} * {dd_ratio}^2 * {ln_dim_time}^2 +"
+           "  {TP[21]} * {dd_ratio}^2 * {number} +"
+           "  {TP[22]} * {dd_ratio}^2 * {number}^2 +"
+           "  {TP[23]} * {dd_ratio}^2 * {ratio} +"
+           "  {TP[24]} * {dd_ratio}^2 * {ratio}^2 +"
+           "  {TP[25]} * {ln_dim_time} * {number} +"
+           "  {TP[26]} * {ln_dim_time} * {number}^2 +"
+           "  {TP[27]} * {ln_dim_time} * {ratio} +"
+           "  {TP[28]} * {ln_dim_time} * {ratio}^2 +"
+           "  {TP[29]} * {ln_dim_time}^2 * {number} +"
+           "  {TP[30]} * {ln_dim_time}^2 * {number}^2 +"
+           "  {TP[31]} * {ln_dim_time}^2 * {ratio} +"
+           "  {TP[32]} * {ln_dim_time}^2 * {ratio}^2 +"
+           "  {TP[33]} * {number} * {ratio} +"
+           "  {TP[34]} * {number} * {ratio}^2 +"
+           "  {TP[35]} * {number}^2 * {ratio} +"
+           "  {TP[36]} * {number}^2 * {ratio}^2))")
+    rcmd = res.format(out=out, yearly_load=yearly_load,
+                      ground_conductivity=ground_conductivity, length=length,
+                      dd_ratio=dd_ratio, ln_dim_time=ln_dim_time,
+                      number=number, ratio=ratio, TP=TP, pi=pi)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def _get_length(yearly, long_term, monthly, medium_term, peak, short_term,
+                borehole_res, ftemp, gtemp, temp_penality):
+    """Return a first attempt of the BHE field lenght
+
+    Example
+    -------
+    >>> _get_length(yearly=-1762, long_term=0.16963475237211000,
+    ...             monthly=-100000, medium_term=0.16000118270019600,
+    ...             peak=-392250, short_term=0.10067477057511400,
+    ...             borehole_res=0.1015426274842840,
+    ...             ftemp=2.750810810810810, gtemp=12.41,
+    ...             temp_penality=-0.2400252902422710)
+    10151.515582310705
+    """
+    # print(dict(yearly=yearly, long_term=long_term, monthly=monthly,
+    #            medium_term=medium_term, peak=peak, short_term=short_term,
+    #            borehole_res=borehole_res, ftemp=ftemp, gtemp=gtemp,
+    #            temp_penality=temp_penality))
+    return ((yearly * long_term +
+             monthly * medium_term +
+             peak * short_term + peak * borehole_res) /
+            (ftemp - gtemp - temp_penality))
+
+
+def _r_get_length(out, yearly_load, long_term, monthly_load, medium_term,
+                  peak_load, short_term, borehole_resistence,
+                  fluid_temp, ground_temp, temperature_penality,
+                  execute=True, **kwargs):
+    """Return a first attempt of the BHE field lenght
+
+    Example
+    -------
+    >>> _r_get_length('length', yearly_load='y_load', long_term='l_term',
+    ...               monthly_load='m_load', medium_term='m_term',
+    ...               peak_load='p_load', short_term='s_term',
+    ...               borehole_resistence='bh_resistence',
+    ...               fluid_temp='f_temp', ground_temp='g_temp',
+    ...               temperature_penality='temp_pen', execute=False)
+    ...                                        # doctest: +NORMALIZE_WHITESPACE
+    'length = ((y_load * l_term +
+                m_load * m_term +
+                p_load * s_term +
+                p_load * bh_resistence) /
+               (f_temp - g_temp - temp_pen))'
+    """
+    res = ("{out} = (({yearly_load} * {long_term} +"
+           "          {monthly_load} * {medium_term} +"
+           "          {peak_load} * {short_term} +"
+           "          {peak_load} * {borehole_resistence}) /"
+           "         ({fluid_temp} - {ground_temp} - {temp_penality}))")
+    rcmd = res.format(out=out, yearly_load=yearly_load, long_term=long_term,
+                      monthly_load=monthly_load, medium_term=medium_term,
+                      peak_load=peak_load, short_term=short_term,
+                      borehole_resistence=borehole_resistence,
+                      fluid_temp=fluid_temp, ground_temp=ground_temp,
+                      temp_penality=temperature_penality)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+    return rcmd
+
+
+def temperature_penality(field, length):
+    """Return the temperature penalty: Tp [°C]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly=-392250,
+    ...                                    monthly=-100000,
+    ...                                    yearly=-1762),
+    ...           ground=GroundProperties(conductivity=2.25,
+    ...                                   diffusivity=0.06757039008,
+    ...                                   temperature=12.41),
+    ...           fluid=FluidProperties(capacity=4000, massflow=0.074,
+    ...                                 inlettemp=4.44),
+    ...           borehole=Borehole(radius=0.054,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...                             convection=1000.)
+    ...           )
+    >>> field = BoreholeField(distance=6.1, number=120, ratio=1.2, bhe=bhe)
+    >>> temperature_penality(field, 9899.25626464766)
+    -0.24002529024227098
+    """
+    dd_ratio = distance_depth_ratio(field, length)
+    ln_dim_time = log_dimless_time(field, length)
+    return _temperature_penality(field.bhe.ground_loads.yearly,
+                                 field.bhe.ground.conductivity,
+                                 length, dd_ratio, ln_dim_time,
+                                 field.number, field.ratio)
+
+
+def r_temperature_penality(out, field, length, execute=True, **kwargs):
+    """Return the temperature penalty: Tp [°C]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly=-392250,
+    ...                                    monthly=-100000,
+    ...                                    yearly=-1762),
+    ...           ground=GroundProperties(conductivity=2.25,
+    ...                                   diffusivity=0.06757039008,
+    ...                                   temperature=12.41),
+    ...           fluid=FluidProperties(capacity=4000, massflow=0.074,
+    ...                                 inlettemp=4.44),
+    ...           borehole=Borehole(radius=0.054,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...                             convection=1000.)
+    ...           )
+    >>> field = BoreholeField(distance=6.1, number=120, ratio=1.2, bhe=bhe)
+    >>> r_temperature_penality('temp_pen', field, 'length', execute=False)
+    ...                                                    # doctest: +ELLIPSIS
+    'temp_pen = (...)'
+    """
+    dd_ratio = BASENAME + 'dd_ratio'
+    r_distance_depth_ratio(dd_ratio, field, length, execute=execute, **kwargs)
+    dimless_time = BASENAME + 'ln_dimless_time'
+    r_log_dimless_time(dimless_time, field, length, execute=execute, **kwargs)
+    return _r_temperature_penality(out,
+                                   field.bhe.ground_loads.yearly,
+                                   field.bhe.ground.conductivity,
+                                   length, dd_ratio, dimless_time,
+                                   field.number, field.ratio,
+                                   execute=execute, **kwargs)
+
+
+def field_length(field, tol=1e-3):
+    """Return the total borefield length: L [m]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly=-392250,
+    ...                                    monthly=-100000,
+    ...                                    yearly=-1762),
+    ...           ground=GroundProperties(conductivity=2.25,
+    ...                                   diffusivity=0.06757039008,
+    ...                                   temperature=12.41),
+    ...           fluid=FluidProperties(capacity=4000, massflow=0.074,
+    ...                                 inlettemp=4.44),
+    ...           borehole=Borehole(radius=0.054,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...                             convection=1000.)
+    ...           )
+    >>> field = BoreholeField(distance=6.1, number=120, ratio=1.2, bhe=bhe)
+    >>> field_length(field)
+    10149.682910753265
+    """
+    # extract common variables
+    long_term = ground_resistence(field.bhe.ground, field.bhe.borehole.radius, period='10y')
+    medium_term = ground_resistence(field.bhe.ground, field.bhe.borehole.radius, period='1m')
+    short_term = ground_resistence(field.bhe.ground, field.bhe.borehole.radius, period='6h')
+    borehole_res = bh_resistence(field.bhe.borehole, field.bhe.ground.conductivity)
+    peak = field.bhe.ground_loads.hourly
+    monthly = field.bhe.ground_loads.monthly
+    yearly = field.bhe.ground_loads.yearly
+    ftemp = fluid_temperature_borehole(field.bhe.fluid, peak)
+    gtemp = field.bhe.ground.temperature
+    # start iteration get length for a single BHE
+    length0 = bhe_length(field.bhe)
+
+    # get correct length considering the bhe interferences
+    length1 = _get_length(yearly, long_term,
+                          monthly, medium_term, peak, short_term,
+                          borehole_res, ftemp, gtemp,
+                          temperature_penality(field, length0))
+    while abs(length1 - length0) > tol:
+        length0 = length1
+        length1 = _get_length(yearly, long_term,
+                              monthly, medium_term, peak, short_term,
+                              borehole_res, ftemp, gtemp,
+                              temperature_penality(field, length0))
+    return length1
+
+
+def abs_diff_gt_tol(out, raster_a, raster_b, tol=1e-3, execute=True, **kwargs):
+    res = ("{out} = if(abs({raster_a} - {raster_b}) > {tol}, 1, 0)")
+    rcmd = res.format(out=out, raster_a=raster_a, raster_b=raster_b, tol=tol)
+    if execute:
+        grast.mapcalc(rcmd, **kwargs)
+        info = grast.raster_info(out)
+    else:
+        info = dict(max=0)
+    return True if info['max'] == 1 else False
+
+
+def r_field_length(out, field, infovars,
+                   basename=BASENAME,
+                   length_single=None, tol=1e-3, execute=True, **kwargs):
+    """Return the total borefield length: L [m]
+
+    Example
+    -------
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly=-392250,
+    ...                                    monthly=-100000,
+    ...                                    yearly=-1762),
+    ...           ground=GroundProperties(conductivity=2.25,
+    ...                                   diffusivity=0.06757039008,
+    ...                                   temperature=12.41),
+    ...           fluid=FluidProperties(capacity=4000, massflow=0.074,
+    ...                                 inlettemp=4.44),
+    ...           borehole=Borehole(radius=0.054,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.45, k_grout=1.73, distance=0.0471,
+    ...                             convection=1000.)
+    ...           )
+    >>> field = BoreholeField(distance=6.1, number=120, ratio=1.2, bhe=bhe)
+    >>> infovars = get_vars('bhe_field_length', bhe, BASENAME, execute=False)
+    >>> r_field_length('bhe_field_length', field, infovars, execute=False)
+    'bhe_field_length'
+    """
+    # extract common variables
+    peak = field.bhe.ground_loads.hourly
+    monthly = field.bhe.ground_loads.monthly
+    yearly = field.bhe.ground_loads.yearly
+    gtemp = field.bhe.ground.temperature
+    # start iteration get length for a single BHE
+    len_template = basename + out + '_lenght_{:02d}'
+    length0 = (len_template.format(0)
+               if length_single is None else length_single)
+    if not exists(length0):
+        r_bhe_length(length0, field.bhe, infovars, execute=execute, **kwargs)
+
+    temp_penality = basename + out + '_temp_penality'
+    r_temperature_penality(temp_penality, field, length0,
+                           execute=execute, **kwargs)
+
+    # get correct length considering the bhe interferences
+    length1 = len_template.format(1)
+    _r_get_length(length1, yearly, infovars.long_term,
+                  monthly, infovars.medium_term,
+                  peak, infovars.short_term,
+                  infovars.resistence, infovars.fluid_temp,
+                  gtemp, temp_penality, execute=execute, **kwargs)
+    diff_template = basename + out + '_absdiff_{:02d}'
+    index = 1
+    diff = diff_template.format(index)
+    while abs_diff_gt_tol(diff, length0, length1, tol=tol,
+                          execute=execute, **kwargs):
+        index += 1
+        length0 = length1
+        length1 = len_template.format(index)
+        _r_get_length(length1, yearly, infovars.long_term,
+                      monthly, infovars.medium_term,
+                      peak, infovars.short_term,
+                      infovars.resistence, infovars.fluid_temp,
+                      gtemp, temp_penality, execute=execute, **kwargs)
+    rename('raster', length1, out)
+    return out
+
+
+if __name__ == '__main__':
+    import doctest
+    doctest.testmod()


Property changes on: grass-addons/grass7/raster/r.green/r.green.gshp/libgshp/ashrae.py
___________________________________________________________________
Added: svn:mime-type
   + text/x-python
Added: svn:eol-style
   + native

Modified: grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.html
===================================================================
--- grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.html	2017-05-22 06:11:21 UTC (rev 71102)
+++ grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.html	2017-05-22 06:18:32 UTC (rev 71103)
@@ -5,28 +5,17 @@
 
 The module r.green.gshp consists of the following different parts:<br><br>
 
-- <a href="r.green.gshp.theoretical.html">r.green.hydro.theoretical</a><br>
-  calculates for each basin the theoretical maximum hydropower energy potential<br>
-  input raster maps:    - discharge along river network      - elevation of the considered region<br>  
-  input vector map:       existing plant position<br>
-  output vector maps: - available river segments               - optimal plant position<br><br>
+- <a href="r.green.gshp.theoretical.html">r.green.gshp.theoretical</a><br>
+  [status: testing] calculates for each single Borehole Heat Exchanger (BHE) the maximum power/energy potential, assessed through the G.POT method.<br><br>
 
-- <a href="r.green.gshp.planning.html">r.green.hydro.recommended</a><br>
-  detects the potential plant position considering legal and ecological constraints and the user's recommendations<br>
-  input raster maps:    - discharge along river network      - elevation of the considered region     - minimum flow discharge<br>  
-  input vector maps:    - existing plant position                    -  areas excluded from calculation<br>
-  output vector maps: - available river segments               - optimal plant position<br><br>
+- <a href="r.green.gshp.planning.html">r.green.gshp.planning</a><br>
+  [status: not implemented yet] consider legal, ecological constraints and the user's recommendations to plan the use of the GSHP system<br><br>
 
-- <a href="r.green.gshp.technical.html">r.green.hydro.technical</a><br>
-  calculates the hydropower potential considering technical constrains (head losses, efficiency of turbine)<br>
-  input vector map:       intakes and restitutions of the potential plants<br>
-  output vector map:  - structure (derivation channel and penstock) on both sides of the river for each potential plant<br><br>
+- <a href="r.green.gshp.technical.html">r.green.gshp.technical</a><br>
+  [status: testing] calculates the length of a BHE following the ASHRAE methodology.<br><br>
 
-- <a href="r.green.gshp.financial.html">r.green.hydro.financial</a><br>
-  computes the economic costs and values of the plants<br>
-  input raster maps:    - landuse      - slope<br>  
-  input vector maps:   - segments of the potential plants   - plant structure     - electric grid<br>
-  output vector map:    structure of the potential plants with a re-organized table with e.g. power, gross head, total cost<br><br>
+- <a href="r.green.gshp.financial.html">r.green.gshp.financial</a><br>
+  [status: not implemented yet] computes the economic and financial feasibility of the ground source heat plants<br><br>
 
 <h2>AUTHORS</h2>
 

Added: grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/Makefile
===================================================================
--- grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/Makefile	                        (rev 0)
+++ grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/Makefile	2017-05-22 06:18:32 UTC (rev 71103)
@@ -0,0 +1,7 @@
+MODULE_TOPDIR = ../../../..
+
+PGM = r.green.gshp.technical
+
+include $(MODULE_TOPDIR)/include/Make/Script.make
+
+default: script

Added: grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.html
===================================================================
--- grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.html	                        (rev 0)
+++ grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.html	2017-05-22 06:18:32 UTC (rev 71103)
@@ -0,0 +1,25 @@
+<h2>DESCRIPTION</h2>
+<em>r.green.gshp.technical</em> calculates the lenght of the borehole for Groud Source Heat Pump plants.<br><br>
+
+<h2>NOTES</h2>
+
+The required inputs are the ground source conductivity and the ground loads.<br><br>
+
+
+<h2>EXAMPLES</h2>
+
+
+<div class="code"><pre>python r.green.gshp.technical.py --overwrite ground_conductivity=g_conductivity ground_diffusivity_rast=g_diffusivity ground_temp_rast=g_temperature g_loads_6h_rast=g_loads_6h g_loads_1m_rast=g_loads_1m g_loads_1y_rast=g_loads_1y fluid_capacity=4000. fluid_massflow=0.074 fluid_inlettemp=4.44 bh_radius=0.054 pipe_inner_radius=0.01365 pipe_outer_radius=0.0167 bh_convection=1000. pipe_distance=0.0471 k_pipe=0.45 k_grout=1.73 field_distance=6.1 field_number=120 field_ratio=1.2 bhe_length=len_bhe bhe_field_length=len_bhe_field -d --o</pre></div><br>
+
+
+<h2>SEE ALSO</h2>
+<em>
+<a href="r.green.gshp.technical.html">r.green.gshp.technical</a><br>
+</em>
+
+<h2>AUTHORS</h2>
+Pietro Zambelli (Eurac Research, Bolzano, Italy)<br>
+
+<h2>REFERENCES</h2>
+
+<p><i>Last changed: $Date: 2017-04-26 14:53:25 +0200 (Wed, 26 Apr 2017) $</i>

Added: grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.py
===================================================================
--- grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.py	                        (rev 0)
+++ grass-addons/grass7/raster/r.green/r.green.gshp/r.green.gshp.technical/r.green.gshp.technical.py	2017-05-22 06:18:32 UTC (rev 71103)
@@ -0,0 +1,443 @@
+#!/usr/bin/env python2
+# -*- coding: utf-8 -*-
+#
+############################################################################
+#
+# MODULE:      r.green.gshp.technical
+# AUTHOR(S):   Pietro Zambelli
+# PURPOSE:     Calculate the ground source heat pump technical poential using
+#              the ASHRAE method
+# COPYRIGHT:   (C) 2017 by the GRASS Development Team
+#
+#              This program is free software under the GNU General Public
+#              License (>=v2). Read the file COPYING that comes with GRASS
+#              for details.
+#
+#############################################################################
+#
+
+#%module
+#% description: Calculate the Ground Source Heat Pump technical potential using the ASHRAE method.
+#% keywords: raster
+#%end
+
+##
+## REQUIRED INPUTS
+##
+#%option G_OPT_R_INPUT
+#% key: ground_conductivity
+#% description: Raster with depth-averaged ground thermal conductivity λ [W m-1 K-1]
+#% required: yes
+#%end
+
+##
+## OPTIONAL INPUTS
+##
+
+#######################
+## Ground properties ##
+#######################
+#%option G_OPT_R_INPUT
+#% key: ground_diffusivity_rast
+#% description: Raster with depth-averaged ground diffusivity [m2 day-1]
+#% required: no
+#% guisection: Ground
+#%end
+#%option
+#% key: ground_diffusivity_value
+#% type: double
+#% key_desc: double
+#% description: Value with depth-averaged ground diffusivity  [m2 day-1]
+#% required: no
+#% answer: 0.086
+#% guisection: Ground
+#%end
+
+#%option G_OPT_R_INPUT
+#% key: ground_temp_rast
+#% description: Raster with the initial ground temperature T0 [°C]
+#% required: no
+#% guisection: Ground
+#%end
+#%option
+#% key: ground_temp_value
+#% type: double
+#% key_desc: double
+#% description: Value with the initial ground temperature T0 [°C]
+#% required: no
+#% answer: 10.
+#% guisection: Ground
+#%end
+
+
+#############################
+## Ground loads properties ##
+#############################
+#%option G_OPT_R_INPUT
+#% key: g_loads_6h_rast
+#% description: Peak of the maximum 6 hourly ground loads [W]
+#% required: no
+#% guisection: Ground loads
+#%end
+#%option
+#% key: g_loads_6h_value
+#% type: double
+#% key_desc: double
+#% description: Peak of the maximum 6 hourly ground loads [W]
+#% required: no
+#% answer: 60.
+#% guisection: Ground loads
+#%end
+
+#%option G_OPT_R_INPUT
+#% key: g_loads_1m_rast
+#% description: Month with the maximum ground loads [W]
+#% required: no
+#% guisection: Ground loads
+#%end
+#%option
+#% key: g_loads_1m_value
+#% type: double
+#% key_desc: double
+#% description: Month with the maximum ground loads [W]
+#% required: no
+#% answer: 15.
+#% guisection: Ground loads
+#%end
+
+#%option G_OPT_R_INPUT
+#% key: g_loads_1y_rast
+#% description: Yearly average ground loads [W]
+#% required: no
+#% guisection: Ground loads
+#%end
+#%option
+#% key: g_loads_1y_value
+#% type: double
+#% key_desc: double
+#% description: Yearly average ground loads [W]
+#% required: no
+#% answer: 5.
+#% guisection: Ground loads
+#%end
+
+######################
+## Fluid properties ##
+######################
+#%option
+#% key: fluid_capacity
+#% type: double
+#% key_desc: double
+#% description: Fluid capacity Cp [J kg-1 K-1]
+#% required: no
+#% answer: 4200.
+#% guisection: Fluid
+#%end
+#%option
+#% key: fluid_massflow
+#% type: double
+#% key_desc: double
+#% description: Fluid massflow  [kg s-1 kW-1]
+#% required: no
+#% answer: 0.050
+#% guisection: Fluid
+#%end
+#%option
+#% key: fluid_inlettemp
+#% type: double
+#% key_desc: double
+#% description: Inlet temperature  [°C]
+#% required: no
+#% answer: 2.
+#% guisection: Fluid
+#%end
+
+
+#########################
+## Borehole properties ##
+#########################
+#%option
+#% key: bh_radius
+#% type: double
+#% key_desc: double
+#% description: Borehole radius [m]
+#% required: no
+#% answer: 0.06
+#% guisection: Borehole
+#%end
+#%option
+#% key: bh_convection
+#% type: double
+#% key_desc: double
+#% description: Internal convection coefficient [W m-2 K-1]
+#% required: no
+#% answer: 1.52
+#% guisection: Borehole
+#%end
+#%option
+#% key: bh_resistence
+#% type: double
+#% key_desc: double
+#% description: Borehole thermal resistence [m K W-1]
+#% required: no
+#% answer: nan
+#% guisection: Borehole
+#%end
+
+#%option
+#% key: pipe_inner_radius
+#% type: double
+#% key_desc: double
+#% description: Borehole pipe inner radius [m]
+#% required: no
+#% answer: 0.01365
+#% guisection: Borehole
+#%end
+#%option
+#% key: pipe_outer_radius
+#% type: double
+#% key_desc: double
+#% description: Borehole pipe outer radius [m]
+#% required: no
+#% answer: 0.0167
+#% guisection: Borehole
+#%end
+#%option
+#% key: pipe_distance
+#% type: double
+#% key_desc: double
+#% description: Center-to-center distance between pipes [m]
+#% required: no
+#% answer: 0.0511
+#% guisection: Borehole
+#%end
+#%option
+#% key: k_pipe
+#% type: double
+#% key_desc: double
+#% description: Pipe thermal conductivity [W m-1 K-1]
+#% required: no
+#% answer: 0.42
+#% guisection: Borehole
+#%end
+#%option
+#% key: k_grout
+#% type: double
+#% key_desc: double
+#% description: Grout thermal conductivity [W m-1 K-1]
+#% required: no
+#% answer: 1.52
+#% guisection: Borehole
+#%end
+
+
+
+###############################
+## Borehole filed properties ##
+###############################
+#%option
+#% key: field_distance
+#% type: double
+#% key_desc: double
+#% description: Distance between boreholes heat exchanger
+#% required: no
+#% answer: 6.1
+#% guisection: BHE Field
+#%end
+#%option
+#% key: field_number
+#% type: integer
+#% key_desc: integer
+#% description: Number of borehole heat exchanger
+#% required: no
+#% answer: 2
+#% guisection: BHE Field
+#%end
+#%option
+#% key: field_ratio
+#% type: double
+#% key_desc: double
+#% description: Borefield aspect ratio
+#% required: no
+#% answer: 1.2
+#% guisection: BHE Field
+#%end
+
+
+#############
+## Outputs ##
+#############
+#%option G_OPT_R_OUTPUT
+#% key: bhe_length
+#% type: string
+#% key_desc: name
+#% description: Name of output raster map with the geothermal length of the BHE [m]
+#% required: no
+#%end
+#%option G_OPT_R_OUTPUT
+#% key: bhe_field_length
+#% type: string
+#% key_desc: name
+#% description: Name of output raster map with the geothermal length of the BHE field [m]
+#% required: no
+#%end
+
+#%flag
+#% key: d
+#% description: Debug with intermediate maps
+#%end
+
+from __future__ import print_function
+
+import os
+import sys
+import atexit
+
+# import grass libraries
+from grass.script import core as gcore
+from grass.pygrass.utils import set_path
+
+try:
+    # set python path to the shared r.green libraries
+    set_path('r.green', 'libgshp', '..')
+    set_path('r.green', 'libgreen', os.path.join('..', '..'))
+    from libgreen.utils import cleanup, rast_or_numb
+    from libgshp.ashrae import (GroundProperties, GroundLoads,
+                                FluidProperties, Borehole,
+                                BoreholeExchanger, BoreholeField,
+                                get_vars, r_bhe_length, r_field_length)
+except ImportError:
+    try:
+        set_path('r.green', 'libgshp', os.path.join('..', 'etc', 'r.green'))
+        set_path('r.green', 'libgreen', os.path.join('..', 'etc', 'r.green'))
+        from libgreen.utils import cleanup, rast_or_numb
+        from libgshp.ashrae import (GroundProperties, GroundLoads,
+                                    FluidProperties, Borehole,
+                                    BoreholeExchanger, BoreholeField,
+                                    get_vars, r_bhe_length, r_field_length)
+    except ImportError:
+        gcore.warning('libgreen and libgshp not in the python path!')
+
+
+def main(opts, flgs):
+    """
+    Parameters
+    ----------
+
+    ground_conductivity: k [W m-1 K-1]
+        thermal conductivity
+    ground_diffusivity: α [m2 day-1]
+        thermal diffusivity
+    ground_temperature: Tg [°C]
+        Undisturbed ground temperature
+    g_loads_6h [W]
+        Peak of 6 hours ground load
+    g_loads_1m:  [W]
+        Maximum monthly ground load
+    g_loads_1y: qy [W]
+        Average ground load
+    fluid_capacity: Cp [J kg-1 K-1]
+        Thermal heat capacity of the fluid
+    fluid_massflow: mfls [kg s-1 kW-1]
+        total mass flow rate per kW of peak hourly ground load
+    fluid_inlettemp: TinHP [°C]
+        max/min heat pump inlet temperature
+    bh_radius: radius [m]
+        borehole radius
+    bh_convection: hconv [W m-2 K-1]
+        Internal convection coefficient
+    pipe_inner: rpin [m]
+        pipe inner radius
+    pipe_outer: rpext [m]
+        pipe outer radius
+    pipe_distance: LU [m]
+        center-to-center distance between pipes
+    k_pipe: kpipe [W m-1 K-1]
+        pipe thermal conductivity
+    k_grout: kgrout [W m-1 K-1]
+        grout thermal conductivity
+    distance_bhe: B [m]
+        distance between boreholes (BHE)
+    number_bhe: [n]
+        number of boreholes (BHE)
+    ratio_bhe: [-]
+        borefield aspect ratio
+
+    >>> bhe = BoreholeExchanger(
+    ...           ground_loads=GroundLoads(hourly='g_loads_6h',
+    ...                                    monthly='g_loads_1m',
+    ...                                    yearly='g_loads_1y'),
+    ...           ground=GroundProperties(conductivity='g_conductivity',
+    ...                                   diffusivity='g_diffusivity',
+    ...                                   temperature='g_temperature'),
+    ...           fluid=FluidProperties(capacity=4200, massflow=0.050,
+    ...                                 inlettemp=40.2),
+    ...           borehole=Borehole(radius=0.06,
+    ...                             pipe_inner_radius=0.01365,
+    ...                             pipe_outer_radius=0.0167,
+    ...                             k_pipe=0.42, k_grout=1.5, distance=0.0511,
+    ...                             convection=1000.)
+    ...           )
+    >>> infovars = InfoVars('l_term', 'm_term', 's_term', 'f_temp', 'res')
+    >>> r_bhe_length('bhe_length', bhe, infovars, execute=False)
+    """
+    DEBUG = flags['d']
+    OVER = gcore.overwrite()
+    tmpbase = "tmprgreen_%i" % os.getpid()
+    atexit.register(cleanup, pattern=(tmpbase + '*'), debug=DEBUG)
+
+    # ================================================
+    # GROUND
+    # get raster or scalar value
+    ground = GroundProperties(opts['ground_conductivity'],
+                              rast_or_numb('ground_diffusivity_rast',
+                                           'ground_diffusivity_value', opts),
+                              rast_or_numb('ground_temp_rast',
+                                           'ground_temp_value', opts))
+
+    # ================================================
+    # GROUND LOADS
+    ground_loads = GroundLoads(
+            rast_or_numb('g_loads_6h_rast', 'g_loads_6h_value', opts),
+            rast_or_numb('g_loads_1m_rast', 'g_loads_1m_value', opts),
+            rast_or_numb('g_loads_1y_rast', 'g_loads_1y_value', opts))
+
+    # ================================================
+    # FLUID
+    fluid = FluidProperties(float(opts['fluid_capacity']),
+                            float(opts['fluid_massflow']),
+                            float(opts['fluid_inlettemp']))
+
+    # ================================================
+    # BOREHOLE
+    borehole = Borehole(radius=float(opts['bh_radius']),
+                        pipe_inner_radius=float(opts['pipe_inner_radius']),
+                        pipe_outer_radius=float(opts['pipe_outer_radius']),
+                        k_pipe=float(opts['k_pipe']),
+                        k_grout=float(opts['k_grout']),
+                        distance=float(opts['pipe_distance']),
+                        convection=float(opts['bh_convection']))
+
+    bhe = BoreholeExchanger(ground_loads, ground, fluid, borehole)
+
+    field = BoreholeField(float(opts['field_distance']),
+                          float(opts['field_number']),
+                          float(opts['field_ratio']),
+                          bhe)
+
+    # ================================================
+    # START COMPUTATIONS
+    infovars = get_vars(opts['bhe_length'], bhe, tmpbase, execute=True,
+                        overwrite=OVER)
+
+    r_bhe_length(opts['bhe_length'], bhe, infovars, execute=True,
+                 overwrite=OVER)
+    r_field_length(opts['bhe_field_length'], field, infovars, basename=tmpbase,
+                   length_single=opts['bhe_length'],
+                   execute=True, overwrite=OVER)
+
+
+if __name__ == "__main__":
+    options, flags = gcore.parser()
+    main(options, flags)
+    sys.exit(0)



More information about the grass-commit mailing list