diff --git a/README.md b/README.md
index c82a7314..a89d5b4d 100644
--- a/README.md
+++ b/README.md
@@ -1,6 +1,5 @@
 PySD
 ====
-[![Build Status](https://travis-ci.com/JamesPHoughton/pysd.svg?branch=master)](https://travis-ci.com/JamesPHoughton/pysd)
 [![Coverage Status](https://coveralls.io/repos/github/JamesPHoughton/pysd/badge.svg?branch=master)](https://coveralls.io/github/JamesPHoughton/pysd?branch=master)
 [![Anaconda-Server Badge](https://anaconda.org/conda-forge/pysd/badges/version.svg)](https://anaconda.org/conda-forge/pysd)
 [![PyPI version](https://badge.fury.io/py/pysd.svg)](https://badge.fury.io/py/pysd)
diff --git a/pysd/_version.py b/pysd/_version.py
index 38cf6dbe..fcfdf383 100644
--- a/pysd/_version.py
+++ b/pysd/_version.py
@@ -1 +1 @@
-__version__ = "1.9.1"
+__version__ = "1.10.0"
diff --git a/pysd/py_backend/builder.py b/pysd/py_backend/builder.py
index 36204229..fdcc46bc 100644
--- a/pysd/py_backend/builder.py
+++ b/pysd/py_backend/builder.py
@@ -588,11 +588,6 @@ def build_element(element, subscript_dict):
     # external objecets via .add method
     py_expr_no_ADD = ["ADD" not in py_expr for py_expr in element["py_expr"]]
 
-    if element["subs"][0]:
-        new_subs = utils.make_merge_list(element["subs"], subscript_dict)
-    else:
-        new_subs = None
-
     if sum(py_expr_no_ADD) > 1 and element["kind"] not in [
         "stateful",
         "external",
@@ -611,7 +606,7 @@ def build_element(element, subscript_dict):
                     terse=False)
                 coords = {
                     new_dim: coords[dim]
-                    for new_dim, dim in zip(new_subs, coords)
+                    for new_dim, dim in zip(element["merge_subs"], coords)
                 }
                 dims = list(coords)
                 Imports.add("utils", "rearrange")
@@ -631,12 +626,12 @@ def build_element(element, subscript_dict):
     element["subs_dec"] = ""
     element["subs_doc"] = "None"
 
-    if new_subs:
+    if element["merge_subs"]:
         # We add the list of the subs to the __doc__ of the function
         # this will give more information to the user and make possible
         # to rewrite subscripted values with model.run(params=X) or
         # model.run(initial_condition=(n,x))
-        element["subs_doc"] = "%s" % new_subs
+        element["subs_doc"] = "%s" % element["merge_subs"]
         if element["kind"] in ["component", "setup",
                                "constant", "component_ext_data"]:
             # the decorator is not always necessary as the objects
@@ -644,7 +639,8 @@ def build_element(element, subscript_dict):
             # dimensions always, we should try to reduce to the
             # maximum when we use it
             # re arrange the python object
-            element["subs_dec"] = "@subs(%s, _subscript_dict)" % new_subs
+            element["subs_dec"] =\
+                "@subs(%s, _subscript_dict)" % element["merge_subs"]
             Imports.add("subs")
 
     indent = 8
@@ -743,7 +739,6 @@ def merge_partial_elements(element_list):
                 # Use 'expr' for Vensim models, and 'eqn' for Xmile
                 # (This makes the Vensim equation prettier.)
                 eqn = element["expr"] if "expr" in element else element["eqn"]
-
                 outs[name] = {
                     "py_name": element["py_name"],
                     "real_name": element["real_name"],
@@ -751,6 +746,8 @@ def merge_partial_elements(element_list):
                     "py_expr": [element["py_expr"]],  # in a list
                     "unit": element["unit"],
                     "subs": [element["subs"]],
+                    "merge_subs": element["merge_subs"]
+                    if "merge_subs" in element else None,
                     "lims": element["lims"],
                     "eqn": [eqn.replace(r"\ ", "")],
                     "kind": element["kind"],
@@ -771,7 +768,7 @@ def merge_partial_elements(element_list):
     return list(outs.values())
 
 
-def add_stock(identifier, expression, initial_condition, subs):
+def add_stock(identifier, expression, initial_condition, subs, merge_subs):
     """
     Creates new model element dictionaries for the model elements associated
     with a stock.
@@ -791,6 +788,10 @@ def add_stock(identifier, expression, initial_condition, subs):
         List of strings of subscript indices that correspond to the
         list of expressions, and collectively define the shape of the output.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects
+
     Returns
     -------
     reference: str
@@ -831,6 +832,7 @@ def add_stock(identifier, expression, initial_condition, subs):
                 "kind": "setup",
                 "py_expr": initial_condition,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "doc": "Provides initial conditions for %s function"
                         % identifier,
                 "unit": "See docs for %s" % identifier,
@@ -848,6 +850,7 @@ def add_stock(identifier, expression, initial_condition, subs):
                 "kind": "component",
                 "doc": "Provides derivative for %s function" % identifier,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "unit": "See docs for %s" % identifier,
                 "lims": "None",
                 "eqn": "None",
@@ -868,6 +871,7 @@ def add_stock(identifier, expression, initial_condition, subs):
             "lims": "None",
             "eqn": "None",
             "subs": "",
+            "merge_subs": None,
             "kind": "stateful",
             "arguments": "",
         }
@@ -876,7 +880,8 @@ def add_stock(identifier, expression, initial_condition, subs):
     return "%s()" % py_name, new_structure
 
 
-def add_delay(identifier, delay_input, delay_time, initial_value, order, subs):
+def add_delay(identifier, delay_input, delay_time, initial_value, order,
+              subs, merge_subs):
     """
     Creates code to instantiate a stateful 'Delay' object,
     and provides reference to that object's output.
@@ -913,6 +918,10 @@ def add_delay(identifier, delay_input, delay_time, initial_value, order, subs):
         List of strings of subscript indices that correspond to the
         list of expressions, and collectively define the shape of the output.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects
+
     Returns
     -------
     reference: str
@@ -953,6 +962,7 @@ def add_delay(identifier, delay_input, delay_time, initial_value, order, subs):
                 # exist
                 "py_expr": initial_value,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "doc": "Provides initial conditions for %s function" \
                         % identifier,
                 "unit": "See docs for %s" % identifier,
@@ -970,6 +980,7 @@ def add_delay(identifier, delay_input, delay_time, initial_value, order, subs):
                 "kind": "component",
                 "doc": "Provides input for %s function" % identifier,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "unit": "See docs for %s" % identifier,
                 "lims": "None",
                 "eqn": "None",
@@ -991,6 +1002,7 @@ def add_delay(identifier, delay_input, delay_time, initial_value, order, subs):
             "lims": "None",
             "eqn": "None",
             "subs": "",
+            "merge_subs": None,
             "kind": "stateful",
             "arguments": "",
         }
@@ -1059,6 +1071,7 @@ def add_delay_f(identifier, delay_input, delay_time, initial_value):
         "lims": "None",
         "eqn": "None",
         "subs": "",
+        "merge_subs": None,
         "kind": "stateful",
         "arguments": "",
     }
@@ -1067,7 +1080,7 @@ def add_delay_f(identifier, delay_input, delay_time, initial_value):
 
 
 def add_n_delay(identifier, delay_input, delay_time, initial_value, order,
-                subs):
+                subs, merge_subs):
     """
     Creates code to instantiate a stateful 'DelayN' object,
     and provides reference to that object's output.
@@ -1104,6 +1117,10 @@ def add_n_delay(identifier, delay_input, delay_time, initial_value, order,
         List of strings of subscript indices that correspond to the
         list of expressions, and collectively define the shape of the output.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects
+
     Returns
     -------
     reference: str
@@ -1144,6 +1161,7 @@ def add_n_delay(identifier, delay_input, delay_time, initial_value, order,
                 # exist
                 "py_expr": initial_value,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "doc": "Provides initial conditions for %s function" \
                         % identifier,
                 "unit": "See docs for %s" % identifier,
@@ -1161,6 +1179,7 @@ def add_n_delay(identifier, delay_input, delay_time, initial_value, order,
                 "kind": "component",
                 "doc": "Provides input for %s function" % identifier,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "unit": "See docs for %s" % identifier,
                 "lims": "None",
                 "eqn": "None",
@@ -1183,6 +1202,7 @@ def add_n_delay(identifier, delay_input, delay_time, initial_value, order,
             "lims": "None",
             "eqn": "None",
             "subs": "",
+            "merge_subs": None,
             "kind": "stateful",
             "arguments": "",
         }
@@ -1192,7 +1212,7 @@ def add_n_delay(identifier, delay_input, delay_time, initial_value, order,
 
 
 def add_sample_if_true(identifier, condition, actual_value, initial_value,
-                       subs):
+                       subs, merge_subs):
     """
     Creates code to instantiate a stateful 'SampleIfTrue' object,
     and provides reference to that object's output.
@@ -1217,6 +1237,10 @@ def add_sample_if_true(identifier, condition, actual_value, initial_value,
         List of strings of subscript indices that correspond to the
         list of expressions, and collectively define the shape of the output.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects
+
     Returns
     -------
     reference: str
@@ -1251,6 +1275,7 @@ def add_sample_if_true(identifier, condition, actual_value, initial_value,
             "kind": "setup",  # not specified in the model file, but must exist
             "py_expr": initial_value,
             "subs": subs,
+            "merge_subs": merge_subs,
             "doc": "Provides initial conditions for %s function" % identifier,
             "unit": "See docs for %s" % identifier,
             "lims": "None",
@@ -1269,6 +1294,7 @@ def add_sample_if_true(identifier, condition, actual_value, initial_value,
         "lims": "None",
         "eqn": "None",
         "subs": "",
+        "merge_subs": None,
         "kind": "stateful",
         "arguments": ""
     })
@@ -1277,7 +1303,7 @@ def add_sample_if_true(identifier, condition, actual_value, initial_value,
 
 
 def add_n_smooth(identifier, smooth_input, smooth_time, initial_value, order,
-                 subs):
+                 subs, merge_subs):
     """
     Constructs stock and flow chains that implement the calculation of
     a smoothing function.
@@ -1310,6 +1336,10 @@ def add_n_smooth(identifier, smooth_input, smooth_time, initial_value, order,
         List of strings of subscript indices that correspond to the
         list of expressions, and collectively define the shape of the output
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects
+
     Returns
     -------
     reference: str
@@ -1352,6 +1382,7 @@ def add_n_smooth(identifier, smooth_input, smooth_time, initial_value, order,
                 # exist
                 "py_expr": initial_value,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "doc": "Provides initial conditions for %s function" % \
                        identifier,
                 "unit": "See docs for %s" % identifier,
@@ -1369,6 +1400,7 @@ def add_n_smooth(identifier, smooth_input, smooth_time, initial_value, order,
                 "kind": "component",
                 "doc": "Provides input for %s function" % identifier,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "unit": "See docs for %s" % identifier,
                 "lims": "None",
                 "eqn": "None",
@@ -1390,6 +1422,7 @@ def add_n_smooth(identifier, smooth_input, smooth_time, initial_value, order,
             "lims": "None",
             "eqn": "None",
             "subs": "",
+            "merge_subs": None,
             "kind": "stateful",
             "arguments": "",
         }
@@ -1398,7 +1431,8 @@ def add_n_smooth(identifier, smooth_input, smooth_time, initial_value, order,
     return "%s()" % py_name, new_structure
 
 
-def add_n_trend(identifier, trend_input, average_time, initial_trend, subs):
+def add_n_trend(identifier, trend_input, average_time, initial_trend,
+                subs, merge_subs):
     """
     Constructs Trend object.
 
@@ -1420,6 +1454,10 @@ def add_n_trend(identifier, trend_input, average_time, initial_trend, subs):
         List of strings of subscript indices that correspond to the
         list of expressions, and collectively define the shape of the output.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects
+
     Returns
     -------
     reference: str
@@ -1459,6 +1497,7 @@ def add_n_trend(identifier, trend_input, average_time, initial_trend, subs):
                 # exist
                 "py_expr": initial_trend,
                 "subs": subs,
+                "merge_subs": merge_subs,
                 "doc": "Provides initial conditions for %s function"
                         % identifier,
                 "unit": "See docs for %s" % identifier,
@@ -1480,6 +1519,7 @@ def add_n_trend(identifier, trend_input, average_time, initial_trend, subs):
             "lims": "None",
             "eqn": "None",
             "subs": "",
+            "merge_subs": None,
             "kind": "stateful",
             "arguments": "",
         }
@@ -1525,6 +1565,7 @@ def add_initial(identifier, value):
         "lims": "None",
         "eqn": "None",
         "subs": "",
+        "merge_subs": None,
         "kind": "stateful",
         "arguments": "",
     }
@@ -1533,7 +1574,7 @@ def add_initial(identifier, value):
 
 
 def add_ext_data(identifier, file_name, tab, time_row_or_col, cell, subs,
-                 subscript_dict, keyword):
+                 subscript_dict, merge_subs, keyword):
     """
     Constructs a external object for handling Vensim's GET XLS DATA and
     GET DIRECT DATA functionality.
@@ -1563,6 +1604,10 @@ def add_ext_data(identifier, file_name, tab, time_row_or_col, cell, subs,
         Dictionary describing the possible dimensions of the stock's
         subscripts.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects.
+
     keyword: str
         Data retrieval method ('interpolate', 'look forward', 'hold backward').
 
@@ -1580,7 +1625,7 @@ def add_ext_data(identifier, file_name, tab, time_row_or_col, cell, subs,
 
     coords = utils.simplify_subscript_input(
         utils.make_coord_dict(subs, subscript_dict, terse=False),
-        subscript_dict, return_full=False)
+        subscript_dict, return_full=False, merge_subs=merge_subs)
     keyword = (
         "'%s'" % keyword.strip(":").lower() if isinstance(keyword, str) else
         keyword)
@@ -1614,6 +1659,7 @@ def add_ext_data(identifier, file_name, tab, time_row_or_col, cell, subs,
         "lims": "None",
         "eqn": "None",
         "subs": subs,
+        "merge_subs": merge_subs,
         "kind": kind,
         "arguments": "",
     }
@@ -1621,7 +1667,8 @@ def add_ext_data(identifier, file_name, tab, time_row_or_col, cell, subs,
     return "%s(time())" % external["py_name"], [external]
 
 
-def add_ext_constant(identifier, file_name, tab, cell, subs, subscript_dict):
+def add_ext_constant(identifier, file_name, tab, cell,
+                     subs, subscript_dict, merge_subs):
     """
     Constructs a external object for handling Vensim's GET XLS CONSTANT and
     GET DIRECT CONSTANT functionality.
@@ -1648,6 +1695,10 @@ def add_ext_constant(identifier, file_name, tab, cell, subs, subscript_dict):
         Dictionary describing the possible dimensions of the stock's
         subscripts.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects.
+
     Returns
     -------
     reference: str
@@ -1662,7 +1713,7 @@ def add_ext_constant(identifier, file_name, tab, cell, subs, subscript_dict):
 
     coords = utils.simplify_subscript_input(
         utils.make_coord_dict(subs, subscript_dict, terse=False),
-        subscript_dict, return_full=False)
+        subscript_dict, return_full=False, merge_subs=merge_subs)
     name = utils.make_python_identifier("_ext_constant_%s" % identifier)[0]
 
     # Check if the object already exists
@@ -1692,6 +1743,7 @@ def add_ext_constant(identifier, file_name, tab, cell, subs, subscript_dict):
         "lims": "None",
         "eqn": "None",
         "subs": subs,
+        "merge_subs": merge_subs,
         "kind": kind,
         "arguments": "",
     }
@@ -1700,7 +1752,7 @@ def add_ext_constant(identifier, file_name, tab, cell, subs, subscript_dict):
 
 
 def add_ext_lookup(identifier, file_name, tab, x_row_or_col, cell,
-                   subs, subscript_dict):
+                   subs, subscript_dict, merge_subs):
     """
     Constructs a external object for handling Vensim's GET XLS LOOKUPS and
     GET DIRECT LOOKUPS functionality.
@@ -1730,6 +1782,10 @@ def add_ext_lookup(identifier, file_name, tab, x_row_or_col, cell,
         Dictionary describing the possible dimensions of the stock's
         subscripts.
 
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects.
+
     Returns
     -------
     reference: str
@@ -1744,7 +1800,7 @@ def add_ext_lookup(identifier, file_name, tab, x_row_or_col, cell,
 
     coords = utils.simplify_subscript_input(
         utils.make_coord_dict(subs, subscript_dict, terse=False),
-        subscript_dict, return_full=False)
+        subscript_dict, return_full=False, merge_subs=merge_subs)
     name = utils.make_python_identifier("_ext_lookup_%s" % identifier)[0]
 
     # Check if the object already exists
@@ -1774,6 +1830,7 @@ def add_ext_lookup(identifier, file_name, tab, x_row_or_col, cell,
         "lims": "None",
         "eqn": "None",
         "subs": subs,
+        "merge_subs": merge_subs,
         "kind": kind,
         "arguments": "x",
     }
@@ -1826,6 +1883,7 @@ def add_macro(macro_name, filename, arg_names, arg_vals):
         "lims": "None",
         "eqn": "None",
         "subs": "",
+        "merge_subs": None,
         "kind": "stateful",
         "arguments": "",
     }
diff --git a/pysd/py_backend/utils.py b/pysd/py_backend/utils.py
index f9dab9ab..dd63866f 100644
--- a/pysd/py_backend/utils.py
+++ b/pysd/py_backend/utils.py
@@ -168,7 +168,7 @@ def make_coord_dict(subs, subscript_dict, terse=True):
     return coordinates
 
 
-def make_merge_list(subs_list, subscript_dict):
+def make_merge_list(subs_list, subscript_dict, element=""):
     """
     This is for assisting when building xrmerge. From a list of subscript
     lists returns the final subscript list after mergin. Necessary when
@@ -183,6 +183,10 @@ def make_merge_list(subs_list, subscript_dict):
     subscript_dict: dict
         The full dictionary of subscript names and values.
 
+    element: str (optional)
+        Element name, if given it will be printed with any error or
+        warning message. Default is "".
+
     Returns
     -------
     dims: list
@@ -196,19 +200,33 @@ def make_merge_list(subs_list, subscript_dict):
 
     """
     coords_set = [set() for i in range(len(subs_list[0]))]
-    for subs in subs_list:
-        coords = make_coord_dict(subs, subscript_dict, terse=False)
-        [coords_set[i].update(coords[dim]) for i, dim in enumerate(coords)]
+    coords_list = [
+        make_coord_dict(subs, subscript_dict, terse=False)
+        for subs in subs_list
+    ]
+
+    # update coords set
+    [[coords_set[i].update(coords[dim]) for i, dim in enumerate(coords)]
+     for coords in coords_list]
 
     dims = [None] * len(coords_set)
-    for i, (coord1, coord2) in enumerate(zip(coords, coords_set)):
-        if set(coords[coord1]) == coord2:
+    # create an array with the name of the subranges for all merging elements
+    dims_list = np.array([list(coords) for coords in coords_list]).transpose()
+    indexes = np.arange(len(dims))
+
+    for i, coord2 in enumerate(coords_set):
+        dims1 = [
+            dim for dim in dims_list[i]
+            if dim is not None and set(subscript_dict[dim]) == coord2
+        ]
+        if dims1:
             # if the given coordinate already matches return it
-            dims[i] = coord1
+            dims[i] = dims1[0]
         else:
             # find a suitable coordinate
+            other_dims = dims_list[indexes != i]
             for name, elements in subscript_dict.items():
-                if coord2 == set(elements) and name not in subs_list[0]:
+                if coord2 == set(elements) and name not in other_dims:
                     dims[i] = name
                     break
 
@@ -217,10 +235,11 @@ def make_merge_list(subs_list, subscript_dict):
                 # dimension that completes it
                 for name, elements in subscript_dict.items():
                     if coord2.issubset(set(elements))\
-                      and name not in subs_list[0]:
+                      and name not in other_dims:
                         dims[i] = name
                         warnings.warn(
-                            "\nDimension given by subscripts:"
+                            element
+                            + "\nDimension given by subscripts:"
                             + "\n\t{}\nis incomplete ".format(coord2)
                             + "using {} instead.".format(name)
                             + "\nSubscript_dict:"
@@ -237,14 +256,17 @@ def make_merge_list(subs_list, subscript_dict):
                         subscript_dict[name + str(j)] = elements
                         dims[i] = name + str(j)
                         warnings.warn(
-                            "\nAdding new subscript range to subscript_dict:\n"
+                            element
+                            + "\nAdding new subscript range to"
+                            + " subscript_dict:\n"
                             + name + str(j) + ": " + ', '.join(elements))
                         break
 
             if not dims[i]:
                 # not able to find the correct dimension
                 raise ValueError(
-                    "\nImpossible to find the dimension that contains:"
+                    element
+                    + "\nImpossible to find the dimension that contains:"
                     + "\n\t{}\nFor subscript_dict:".format(coord2)
                     + "\n\t{}".format(subscript_dict)
                 )
@@ -714,8 +736,7 @@ def round_(x):
     return round(x)
 
 
-def simplify_subscript_input(coords, subscript_dict,
-                             return_full=True):
+def simplify_subscript_input(coords, subscript_dict, return_full, merge_subs):
     """
     Parameters
     ----------
@@ -725,13 +746,13 @@ def simplify_subscript_input(coords, subscript_dict,
     subscript_dict: dict
         The subscript dictionary of the model file.
 
-    return_full: bool (optional)
+    return_full: bool
         If True the when coords == subscript_dict, '_subscript_dict'
-        will be returned. Default is True
+        will be returned
 
-    partial: bool (optional)
-        If True "_subscript_dict" will not be returned as possible dict.
-        Used when subscript_dict is not the full dictionary. Default is False.
+    merge_subs: list of strings
+        List of the final subscript range of the python array after
+        merging with other objects
 
     Returns
     -------
@@ -745,14 +766,14 @@ def simplify_subscript_input(coords, subscript_dict,
         return "_subscript_dict"
 
     coordsp = []
-    for dim, coord in coords.items():
+    for ndim, (dim, coord) in zip(merge_subs, coords.items()):
         # find dimensions can be retrieved from _subscript_dict
         if coord == subscript_dict[dim]:
             # use _subscript_dict
-            coordsp.append(f"'{dim}': _subscript_dict['{dim}']")
+            coordsp.append(f"'{ndim}': _subscript_dict['{dim}']")
         else:
             # write whole dict
-            coordsp.append(f"'{dim}': {coord}")
+            coordsp.append(f"'{ndim}': {coord}")
 
     return "{" + ", ".join(coordsp) + "}"
 
diff --git a/pysd/py_backend/vensim/vensim2py.py b/pysd/py_backend/vensim/vensim2py.py
index e9f0eda6..f3c5bdb3 100644
--- a/pysd/py_backend/vensim/vensim2py.py
+++ b/pysd/py_backend/vensim/vensim2py.py
@@ -565,7 +565,7 @@ def visit_var_definition(self, n, vc):
 
         def generic_visit(self, n, vc):
             return "".join(filter(None, vc)) or n.text or ""
-        
+
     tree = parser.parse(sketch_line)
     return SketchParser(tree, namespace=namespace).view_or_var
 
@@ -774,114 +774,141 @@ def parse_units(units_str):
 }
 
 builders = {
-    "integ": lambda element, subscript_dict, args: builder.add_stock(
+    "integ": lambda element, subscript_dict, merge_subs, args:
+    builder.add_stock(
         identifier=element["py_name"],
         expression=args[0],
         initial_condition=args[1],
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "delay1": lambda element, subscript_dict, args: builder.add_delay(
+    "delay1": lambda element, subscript_dict, merge_subs, args:
+    builder.add_delay(
         identifier=element["py_name"],
         delay_input=args[0],
         delay_time=args[1],
         initial_value=args[0],
         order="1",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "delay1i": lambda element, subscript_dict, args: builder.add_delay(
+    "delay1i": lambda element, subscript_dict, merge_subs, args:
+    builder.add_delay(
         identifier=element["py_name"],
         delay_input=args[0],
         delay_time=args[1],
         initial_value=args[2],
         order="1",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "delay3": lambda element, subscript_dict, args: builder.add_delay(
+    "delay3": lambda element, subscript_dict, merge_subs, args:
+    builder.add_delay(
         identifier=element["py_name"],
         delay_input=args[0],
         delay_time=args[1],
         initial_value=args[0],
         order="3",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "delay3i": lambda element, subscript_dict, args: builder.add_delay(
+    "delay3i": lambda element, subscript_dict, merge_subs, args:
+    builder.add_delay(
         identifier=element["py_name"],
         delay_input=args[0],
         delay_time=args[1],
         initial_value=args[2],
         order="3",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "delay fixed": lambda element, subscript_dict, args: builder.add_delay_f(
+    "delay fixed": lambda element, subscript_dict, merge_subs, args:
+    builder.add_delay_f(
         identifier=element["py_name"],
         delay_input=args[0],
         delay_time=args[1],
-        initial_value=args[2],
+        initial_value=args[2]
     ),
-    "delay n": lambda element, subscript_dict, args: builder.add_n_delay(
+    "delay n": lambda element, subscript_dict, merge_subs, args:
+    builder.add_n_delay(
         identifier=element["py_name"],
         delay_input=args[0],
         delay_time=args[1],
         initial_value=args[2],
         order=args[3],
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "sample if true": lambda element, subscript_dict, args:
+    "sample if true": lambda element, subscript_dict, merge_subs, args:
     builder.add_sample_if_true(
         identifier=element["py_name"],
         condition=args[0],
         actual_value=args[1],
         initial_value=args[2],
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "smooth": lambda element, subscript_dict, args: builder.add_n_smooth(
+    "smooth": lambda element, subscript_dict, merge_subs, args:
+    builder.add_n_smooth(
         identifier=element["py_name"],
         smooth_input=args[0],
         smooth_time=args[1],
         initial_value=args[0],
         order="1",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "smoothi": lambda element, subscript_dict, args: builder.add_n_smooth(
+    "smoothi": lambda element, subscript_dict, merge_subs, args:
+    builder.add_n_smooth(
         identifier=element["py_name"],
         smooth_input=args[0],
         smooth_time=args[1],
         initial_value=args[2],
         order="1",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "smooth3": lambda element, subscript_dict, args: builder.add_n_smooth(
+    "smooth3": lambda element, subscript_dict, merge_subs, args:
+    builder.add_n_smooth(
         identifier=element["py_name"],
         smooth_input=args[0],
         smooth_time=args[1],
         initial_value=args[0],
         order="3",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "smooth3i": lambda element, subscript_dict, args: builder.add_n_smooth(
+    "smooth3i": lambda element, subscript_dict, merge_subs, args:
+    builder.add_n_smooth(
         identifier=element["py_name"],
         smooth_input=args[0],
         smooth_time=args[1],
         initial_value=args[2],
         order="3",
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "smooth n": lambda element, subscript_dict, args: builder.add_n_smooth(
+    "smooth n": lambda element, subscript_dict, merge_subs, args:
+    builder.add_n_smooth(
         identifier=element["py_name"],
         smooth_input=args[0],
         smooth_time=args[1],
         initial_value=args[2],
         order=args[3],
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "trend": lambda element, subscript_dict, args: builder.add_n_trend(
+    "trend": lambda element, subscript_dict, merge_subs, args:
+    builder.add_n_trend(
         identifier=element["py_name"],
         trend_input=args[0],
         average_time=args[1],
         initial_trend=args[2],
         subs=element["subs"],
+        merge_subs=merge_subs
     ),
-    "get xls data": lambda element, subscript_dict, args: builder.add_ext_data(
+    "get xls data": lambda element, subscript_dict, merge_subs, args:
+    builder.add_ext_data(
         identifier=element["py_name"],
         file_name=args[0],
         tab=args[1],
@@ -889,9 +916,10 @@ def parse_units(units_str):
         cell=args[3],
         subs=element["subs"],
         subscript_dict=subscript_dict,
+        merge_subs=merge_subs,
         keyword=element["keyword"],
     ),
-    "get xls constants": lambda element, subscript_dict, args:
+    "get xls constants": lambda element, subscript_dict, merge_subs, args:
     builder.add_ext_constant(
         identifier=element["py_name"],
         file_name=args[0],
@@ -899,8 +927,9 @@ def parse_units(units_str):
         cell=args[2],
         subs=element["subs"],
         subscript_dict=subscript_dict,
+        merge_subs=merge_subs,
     ),
-    "get xls lookups": lambda element, subscript_dict, args:
+    "get xls lookups": lambda element, subscript_dict, merge_subs, args:
     builder.add_ext_lookup(
         identifier=element["py_name"],
         file_name=args[0],
@@ -909,12 +938,13 @@ def parse_units(units_str):
         cell=args[3],
         subs=element["subs"],
         subscript_dict=subscript_dict,
+        merge_subs=merge_subs,
     ),
-    "initial": lambda element, subscript_dict, args:
+    "initial": lambda element, subscript_dict, merge_subs, args:
     builder.add_initial(
         identifier=element["py_name"],
         value=args[0]),
-    "a function of": lambda element, subscript_dict, args:
+    "a function of": lambda element, subscript_dict, merge_subs, args:
     builder.add_incomplete(
         element["real_name"], args
     ),
@@ -1269,16 +1299,14 @@ def visit_array(self, n, vc):
                 else:
                     datastr = n.text
 
-                # Create a cleaner dictionary of coords using _subscript_dict
-                # when possible
-                utils.simplify_subscript_input(coords, subscript_dict,
-                                               return_full=True)
-
                 return builder.build_function_call(
                     functions_utils["DataArray"],
                     [datastr,
-                     utils.simplify_subscript_input(coords, subscript_dict),
-                     repr(list(coords))]
+                     utils.simplify_subscript_input(
+                         coords, subscript_dict,
+                         return_full=True,
+                         merge_subs=element["merge_subs"]),
+                     repr(element["merge_subs"])]
                     )
             else:
                 return n.text.replace(" ", "")
@@ -1339,14 +1367,21 @@ def visit_subscript_list(self, n, vc):
         def visit_build_call(self, n, vc):
             # use only the dict with the final subscripts
             # needed for the good working of externals
-
-            subs = {
+            subs_dict = {
                 k: subscript_dict[k] for k in
-                elements_subs_dict[element["py_name"]]
+                element["merge_subs"]
             }
+            # add subscript ranges given in expr
+            subs_dict.update({
+                sub: subscript_dict[sub] for sub in element['subs']
+                if sub in subscript_dict
+            })
+
             self.kind = "component"
             builder_name = vc[0].strip().lower()
-            name, structure = builders[builder_name](element, subs, vc[4])
+            name, structure = builders[builder_name](
+                element, subs_dict, element["merge_subs"],
+                vc[4])
 
             self.new_structure += structure
 
@@ -1454,8 +1489,20 @@ def visit_regularLookup(self, n, vc):
         def visit_excelLookup(self, n, vc):
             arglist = vc[3].split(",")
             arglist = [arg.replace("\\ ", "") for arg in arglist]
+            # use only the dict with the final subscripts
+            # needed for the good working of externals
+            subs_dict = {
+                k: subscript_dict[k] for k in
+                element["merge_subs"]
+            }
+            # add subscript ranges given in expr
+            subs_dict.update({
+                sub: subscript_dict[sub] for sub in element['subs']
+                if sub in subscript_dict
+            })
             trans, structure = builders["get xls lookups"](
-                element, subscript_dict, arglist
+                element, subs_dict,
+                element["merge_subs"], arglist
             )
 
             self.translation = trans
@@ -1536,10 +1583,17 @@ def translate_section(section, macro_list, sketch, root_path, subview_sep=""):
                 elements_subs_dict[element["py_name"]] = [element["subs"]]
 
     elements_subs_dict = {
-        el: utils.make_merge_list(elements_subs_dict[el], subscript_dict)
+        el: utils.make_merge_list(elements_subs_dict[el], subscript_dict, el)
         for el in elements_subs_dict
     }
 
+    for element in model_elements:
+        if "py_name" in element and element["py_name"] in elements_subs_dict:
+            element["merge_subs"] =\
+                elements_subs_dict[element["py_name"]]
+        else:
+            element["merge_subs"] = None
+
     # Parse components to python syntax.
     for element in model_elements:
         if (element["kind"] == "component" and "py_expr" not in element) or \
@@ -1551,15 +1605,16 @@ def translate_section(section, macro_list, sketch, root_path, subview_sep=""):
                 namespace=namespace,
                 subscript_dict=subscript_dict,
                 macro_list=macro_list,
-                elements_subs_dict=elements_subs_dict,
                 subs_compatibility=subs_compatibility_dict,
+                elements_subs_dict=elements_subs_dict
             )
             element.update(translation)
             model_elements += new_structure
 
         elif element["kind"] == "lookup":
             translation, new_structure = parse_lookup_expression(
-                element, subscript_dict
+                element,
+                subscript_dict=subscript_dict
             )
             element.update(translation)
             model_elements += new_structure
diff --git a/pysd/py_backend/xmile/SMILE2Py.py b/pysd/py_backend/xmile/SMILE2Py.py
index c177fe2f..a1b4d5f8 100644
--- a/pysd/py_backend/xmile/SMILE2Py.py
+++ b/pysd/py_backend/xmile/SMILE2Py.py
@@ -153,7 +153,8 @@
         delay_time=args[1],
         initial_value=args[2] if len(args) > 2 else args[0],
         order="1",
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=None
         ),
 
     "delay3": lambda element, subscript_dict, args: builder.add_n_delay(
@@ -162,7 +163,8 @@
         delay_time=args[1],
         initial_value=args[2] if len(args) > 2 else args[0],
         order="3",
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=None
         ),
 
     "delayn": lambda element, subscript_dict, args: builder.add_n_delay(
@@ -171,7 +173,8 @@
         delay_time=args[1],
         initial_value=args[2] if len(args) > 3 else args[0],
         order=args[2],
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=None
         ),
 
     "smth1": lambda element, subscript_dict, args: builder.add_n_smooth(
@@ -180,7 +183,8 @@
         smooth_time=args[1],
         initial_value=args[2] if len(args) > 2 else args[0],
         order="1",
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=None
         ),
 
     "smth3": lambda element, subscript_dict, args: builder.add_n_smooth(
@@ -189,7 +193,8 @@
         smooth_time=args[1],
         initial_value=args[2] if len(args) > 2 else args[0],
         order="3",
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=None
         ),
 
     "smthn": lambda element, subscript_dict, args: builder.add_n_smooth(
@@ -198,7 +203,8 @@
         smooth_time=args[1],
         initial_value=args[2] if len(args) > 3 else args[0],
         order=args[2],
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=None
         ),
 
     # "forcst" !TODO!
@@ -208,7 +214,8 @@
         trend_input=args[0],
         average_time=args[1],
         initial_trend=args[2] if len(args) > 2 else 0,
-        subs=element["subs"]
+        subs=element["subs"],
+        merge_subs=None
         ),
 
     "init": lambda element, subscript_dict, args: builder.add_initial(
diff --git a/pysd/py_backend/xmile/xmile2py.py b/pysd/py_backend/xmile/xmile2py.py
index b7614d69..7b030144 100644
--- a/pysd/py_backend/xmile/xmile2py.py
+++ b/pysd/py_backend/xmile/xmile2py.py
@@ -248,6 +248,7 @@ def parse_lookup_xml_node(node):
 
         py_expr, new_structure = builder.add_stock(identifier=py_name,
                                                    subs=[],  # Todo later
+                                                   merge_subs=[],
                                                    expression=py_ddt,
                                                    initial_condition=py_initial_value
                                                    )
diff --git a/tests/integration_test_vensim_pathway.py b/tests/integration_test_vensim_pathway.py
index c1d50f38..f044f812 100644
--- a/tests/integration_test_vensim_pathway.py
+++ b/tests/integration_test_vensim_pathway.py
@@ -99,6 +99,13 @@ def test_game(self):
         output, canon = runner(test_models + '/game/test_game.mdl')
         assert_frames_close(output, canon, rtol=rtol)
 
+    def test_get_constants_subrange(self):
+        output, canon = runner(
+          test_models + '/get_constants_subranges/'
+          + 'test_get_constants_subranges.mdl'
+        )
+        assert_frames_close(output, canon, rtol=rtol)
+
     def test_get_data_args_3d_xls(self):
         """
         Test for usage of GET DIRECT/XLS DATA with arguments from a Excel file
diff --git a/tests/test-models b/tests/test-models
index 6b8e96dd..8f9151e4 160000
--- a/tests/test-models
+++ b/tests/test-models
@@ -1 +1 @@
-Subproject commit 6b8e96dd1a4c872e102d98398207abb9293ea72d
+Subproject commit 8f9151e4aefa18cdf51b3c05a3249a47fca04606
diff --git a/tests/unit_test_builder.py b/tests/unit_test_builder.py
index f9c245dc..c4a33d71 100644
--- a/tests/unit_test_builder.py
+++ b/tests/unit_test_builder.py
@@ -23,6 +23,7 @@ def test_no_subs_constant(self):
         string = textwrap.dedent(
             build_element(element={'kind': 'constant',
                                    'subs': [[]],
+                                   'merge_subs': [],
                                    'doc': '',
                                    'py_name': 'my_variable',
                                    'real_name': 'My Variable',
@@ -43,6 +44,7 @@ def test_no_subs_call(self):
         string = textwrap.dedent(
             build_element(element={'kind': 'constant',
                                    'subs': [[]],
+                                   'merge_subs': [],
                                    'doc': '',
                                    'py_name': 'my_first_variable',
                                    'real_name': 'My Variable',
@@ -126,18 +128,30 @@ def test_single_set(self):
 
         self.assertEqual(
             merge_partial_elements(
-                [{'py_name': 'a', 'py_expr': 'ms', 'subs': ['Name1', 'element1'],
-                  'real_name': 'A', 'doc': 'Test', 'unit': None, 'eqn': 'eq1', 'lims': '',
+                [{'py_name': 'a', 'py_expr': 'ms',
+                  'subs': ['Name1', 'element1'],
+                  'merge_subs': ['Name1', 'Elements'],
+                  'real_name': 'A', 'doc': 'Test', 'unit': None,
+                  'eqn': 'eq1', 'lims': '',
                   'kind': 'component', 'arguments': ''},
-                 {'py_name': 'a', 'py_expr': 'njk', 'subs': ['Name1', 'element2'],
-                  'real_name': 'A', 'doc': None, 'unit': None, 'eqn': 'eq2', 'lims': '',
+                 {'py_name': 'a', 'py_expr': 'njk',
+                  'subs': ['Name1', 'element2'],
+                  'merge_subs': ['Name1', 'Elements'],
+                  'real_name': 'A', 'doc': None, 'unit': None,
+                  'eqn': 'eq2', 'lims': '',
                   'kind': 'component', 'arguments': ''},
-                 {'py_name': 'a', 'py_expr': 'as', 'subs': ['Name1', 'element3'],
-                  'real_name': 'A', 'doc': '', 'unit': None, 'eqn': 'eq3', 'lims': '',
+                 {'py_name': 'a', 'py_expr': 'as',
+                  'subs': ['Name1', 'element3'],
+                  'merge_subs': ['Name1', 'Elements'],
+                  'real_name': 'A', 'doc': '', 'unit': None,
+                  'eqn': 'eq3', 'lims': '',
                   'kind': 'component', 'arguments': ''}]),
             [{'py_name': 'a',
               'py_expr': ['ms', 'njk', 'as'],
-              'subs': [['Name1', 'element1'], ['Name1', 'element2'], ['Name1', 'element3']],
+              'subs': [['Name1', 'element1'],
+                       ['Name1', 'element2'],
+                       ['Name1', 'element3']],
+              'merge_subs': ['Name1', 'Elements'],
               'kind': 'component',
               'doc': 'Test',
               'real_name': 'A',
@@ -151,27 +165,36 @@ def test_multiple_sets(self):
         from pysd.py_backend.builder import merge_partial_elements
         actual = merge_partial_elements(
             [{'py_name': 'a', 'py_expr': 'ms', 'subs': ['Name1', 'element1'],
-              'real_name': 'A', 'doc': 'Test', 'unit': None, 'eqn': 'eq1', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'A', 'doc': 'Test', 'unit': None,
+              'eqn': 'eq1', 'lims': '', 'kind': 'component', 'arguments': ''},
              {'py_name': 'a', 'py_expr': 'njk', 'subs': ['Name1', 'element2'],
-              'real_name': 'A', 'doc': None, 'unit': None, 'eqn': 'eq2', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'A', 'doc': None, 'unit': None,
+              'eqn': 'eq2', 'lims': '', 'kind': 'component', 'arguments': ''},
              {'py_name': 'a', 'py_expr': 'as', 'subs': ['Name1', 'element3'],
-              'real_name': 'A', 'doc': '', 'unit': None, 'eqn': 'eq3', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'A', 'doc': '', 'unit': None,
+              'eqn': 'eq3', 'lims': '', 'kind': 'component', 'arguments': ''},
              {'py_name': 'b', 'py_expr': 'bgf', 'subs': ['Name1', 'element1'],
-              'real_name': 'B', 'doc': 'Test', 'unit': None, 'eqn': 'eq4', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'B', 'doc': 'Test', 'unit': None,
+              'eqn': 'eq4', 'lims': '', 'kind': 'component', 'arguments': ''},
              {'py_name': 'b', 'py_expr': 'r4', 'subs': ['Name1', 'element2'],
-              'real_name': 'B', 'doc': None, 'unit': None, 'eqn': 'eq5', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'B', 'doc': None, 'unit': None,
+              'eqn': 'eq5', 'lims': '', 'kind': 'component', 'arguments': ''},
              {'py_name': 'b', 'py_expr': 'ymt', 'subs': ['Name1', 'element3'],
-              'real_name': 'B', 'doc': '', 'unit': None, 'eqn': 'eq6', 'lims': '',
-              'kind': 'component', 'arguments': ''}])
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'B', 'doc': '', 'unit': None,
+              'eqn': 'eq6', 'lims': '', 'kind': 'component', 'arguments': ''}])
 
         expected = [{'py_name': 'a',
                      'py_expr': ['ms', 'njk', 'as'],
-                     'subs': [['Name1', 'element1'], ['Name1', 'element2'], ['Name1', 'element3']],
+                     'subs': [['Name1', 'element1'],
+                              ['Name1', 'element2'],
+                              ['Name1', 'element3']],
+                     'merge_subs': ['Name1', 'Elements'],
                      'kind': 'component',
                      'doc': 'Test',
                      'real_name': 'A',
@@ -182,7 +205,10 @@ def test_multiple_sets(self):
                      },
                     {'py_name': 'b',
                      'py_expr': ['bgf', 'r4', 'ymt'],
-                     'subs': [['Name1', 'element1'], ['Name1', 'element2'], ['Name1', 'element3']],
+                     'subs': [['Name1', 'element1'],
+                              ['Name1', 'element2'],
+                              ['Name1', 'element3']],
+                     'merge_subs': ['Name1', 'Elements'],
                      'kind': 'component',
                      'doc': 'Test',
                      'real_name': 'B',
@@ -198,19 +224,23 @@ def test_non_set(self):
         from pysd.py_backend.builder import merge_partial_elements
         actual = merge_partial_elements(
             [{'py_name': 'a', 'py_expr': 'ms', 'subs': ['Name1', 'element1'],
-              'real_name': 'A', 'doc': 'Test', 'unit': None, 'eqn': 'eq1', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'A', 'doc': 'Test', 'unit': None,
+              'eqn': 'eq1', 'lims': '', 'kind': 'component', 'arguments': ''},
              {'py_name': 'a', 'py_expr': 'njk', 'subs': ['Name1', 'element2'],
-              'real_name': 'A', 'doc': None, 'unit': None, 'eqn': 'eq2', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'Elements'],
+              'real_name': 'A', 'doc': None, 'unit': None,
+              'eqn': 'eq2', 'lims': '', 'kind': 'component', 'arguments': ''},
              {'py_name': 'c', 'py_expr': 'as', 'subs': ['Name1', 'element3'],
-              'real_name': 'C', 'doc': 'hi', 'unit': None, 'eqn': 'eq3', 'lims': '',
-              'kind': 'component', 'arguments': ''},
+              'merge_subs': ['Name1', 'elements3'],
+              'real_name': 'C', 'doc': 'hi', 'unit': None,
+              'eqn': 'eq3', 'lims': '', 'kind': 'component', 'arguments': ''},
              ])
 
         expected = [{'py_name': 'a',
                      'py_expr': ['ms', 'njk'],
                      'subs': [['Name1', 'element1'], ['Name1', 'element2']],
+                     'merge_subs': ['Name1', 'Elements'],
                      'kind': 'component',
                      'doc': 'Test',
                      'real_name': 'A',
@@ -222,6 +252,7 @@ def test_non_set(self):
                     {'py_name': 'c',
                      'py_expr': ['as'],
                      'subs': [['Name1', 'element3']],
+                     'merge_subs': ['Name1', 'elements3'],
                      'kind': 'component',
                      'doc': 'hi',
                      'real_name': 'C',
diff --git a/tests/unit_test_utils.py b/tests/unit_test_utils.py
index 755085ae..bf25ccd3 100644
--- a/tests/unit_test_utils.py
+++ b/tests/unit_test_utils.py
@@ -358,7 +358,7 @@ def test_make_merge_list(self):
             make_merge_list([["dim1", "A", "dim"],
                             ["dim1", "B", "dim"],
                             ["dim1", "C", "dim"]],
-                            subscript_dict),
+                            subscript_dict)
             # use only user warnings
             wu = [w for w in ws if issubclass(w.category, UserWarning)]
             self.assertTrue(len(wu), 1)
@@ -366,6 +366,30 @@ def test_make_merge_list(self):
                 "Adding new subscript range to subscript_dict:\ndim2: A, B, C",
                 str(wu[0].message))
 
+        subscript_dict2 = {
+            "dim1": ["A", "B", "C", "D"],
+            "dim1n": ["A", "B"],
+            "dim1y": ["C", "D"],
+            "dim2": ["E", "F", "G", "H"],
+            "dim2n": ["E", "F"],
+            "dim2y": ["G", "H"]
+        }
+
+        # merging two subranges
+        self.assertEqual(
+            make_merge_list([["dim1y"],
+                             ["dim1n"]],
+                            subscript_dict2),
+            ["dim1"])
+
+        # final subscript in list
+        self.assertEqual(
+            make_merge_list([["dim1", "dim2n"],
+                             ["dim1n", "dim2y"],
+                             ["dim1y", "dim2y"]],
+                            subscript_dict2),
+            ["dim1", "dim2"])
+
     def test_compute_shape(self):
         """"
         Test for computing the shape of an array giving coordinates dictionary
diff --git a/tests/unit_test_vensim2py.py b/tests/unit_test_vensim2py.py
index e15c347a..82bff28f 100644
--- a/tests/unit_test_vensim2py.py
+++ b/tests/unit_test_vensim2py.py
@@ -304,8 +304,9 @@ def test_get_lookup(self):
                 + r"'SheetName', 'index'\ , 'values'))",
                 "py_name": "get_lookup",
                 "subs": [],
+                "merge_subs": []
             },
-            {},
+            {}
         )[1][0]
 
         self.assertEqual(
@@ -477,10 +478,10 @@ def test_stock_construction_function_no_subscripts(self):
             {
                 "expr": "INTEG (FlowA, -10)",
                 "py_name": "test_stock",
-                "subs": []
+                "subs": [],
+                "merge_subs": []
             },
-            {"FlowA": "flowa"},
-            elements_subs_dict={"test_stock": []},
+            {"FlowA": "flowa"}
         )
 
         self.assertEqual(res[1][0]["kind"], "stateful")
@@ -499,13 +500,13 @@ def test_delay_construction_function_no_subscripts(self):
                 "expr": "DELAY1(Variable, DelayTime)",
                 "py_name": "test_delay",
                 "subs": [],
+                "merge_subs": []
             },
             {
                 "Variable": "variable",
                 "DelayTime": "delaytime",
                 "TIME STEP": "time_step",
-            },
-            elements_subs_dict={"test_delay": {}},
+            }
         )
 
         def time_step():
@@ -532,9 +533,9 @@ def test_smooth_construction_function_no_subscripts(self):
                 "expr": "SMOOTH(Variable, DelayTime)",
                 "py_name": "test_smooth",
                 "subs": [],
+                "merge_subs": []
             },
             {"Variable": "variable", "DelayTime": "delaytime"},
-            elements_subs_dict={"test_smooth": []},
         )
 
         # check stateful object creation
@@ -548,9 +549,17 @@ def test_smooth_construction_function_no_subscripts(self):
     def test_subscript_float_initialization(self):
         from pysd.py_backend.vensim.vensim2py import parse_general_expression
 
-        _subscript_dict = {"Dim1": ["A", "B", "C"], "Dim2": ["D", "E"]}
+        _subscript_dict = {
+            "Dim": ["A", "B", "C", "D", "E"],
+            "Dim1": ["A", "B", "C"], "Dim2": ["D", "E"]
+        }
+
+        # case 1
         element = parse_general_expression(
-            {"expr": "3.32", "subs": ["Dim1"]}, {}, _subscript_dict
+            {"expr": "3.32", "subs": ["Dim1"], "py_name": "var",
+             "merge_subs": ["Dim1"]}, {},
+            _subscript_dict
+
         )
         string = element[0]["py_expr"]
         # TODO we should use a = eval(string)
@@ -568,12 +577,35 @@ def test_subscript_float_initialization(self):
         )
         self.assertEqual(a.loc[{"Dim1": "B"}], 3.32)
 
+        # case 2: xarray subscript is a subrange from the final subscript range
+        element = parse_general_expression(
+            {"expr": "3.32", "subs": ["Dim1"], "py_name": "var",
+             "merge_subs": ["Dim"]}, {}, _subscript_dict
+        )
+        string = element[0]["py_expr"]
+        # TODO we should use a = eval(string)
+        # hoewever eval is not detecting _subscript_dict variable
+        self.assertEqual(
+            string,
+            "xr.DataArray(3.32,{'Dim': _subscript_dict['Dim1']},['Dim'])",
+        )
+        a = xr.DataArray(
+            3.32, {"Dim": _subscript_dict["Dim1"]}, ["Dim"]
+        )
+        self.assertDictEqual(
+            {key: list(val.values) for key, val in a.coords.items()},
+            {"Dim": ["A", "B", "C"]},
+        )
+        self.assertEqual(a.loc[{"Dim": "B"}], 3.32)
+
     def test_subscript_1d_constant(self):
         from pysd.py_backend.vensim.vensim2py import parse_general_expression
 
         _subscript_dict = {"Dim1": ["A", "B", "C"], "Dim2": ["D", "E"]}
         element = parse_general_expression(
-            {"expr": "1, 2, 3", "subs": ["Dim1"]}, {}, _subscript_dict
+            {"expr": "1, 2, 3", "subs": ["Dim1"], "py_name": "var",
+             "merge_subs": ["Dim1"]},
+            {}, _subscript_dict
         )
         string = element[0]["py_expr"]
         # TODO we should use a = eval(string)
@@ -597,8 +629,9 @@ def test_subscript_2d_constant(self):
 
         _subscript_dict = {"Dim1": ["A", "B", "C"], "Dim2": ["D", "E"]}
         element = parse_general_expression(
-            {"expr": "1, 2; 3, 4; 5, 6;", "subs": ["Dim1", "Dim2"]}, {},
-            _subscript_dict
+            {"expr": "1, 2; 3, 4; 5, 6;", "subs": ["Dim1", "Dim2"],
+             "merge_subs": ["Dim1", "Dim2"], "py_name": "var"},
+            {}, _subscript_dict
         )
         string = element[0]["py_expr"]
         a = eval(string)
@@ -614,8 +647,9 @@ def test_subscript_3d_depth(self):
 
         _subscript_dict = {"Dim1": ["A", "B", "C"], "Dim2": ["D", "E"]}
         element = parse_general_expression(
-            {"expr": "1, 2; 3, 4; 5, 6;", "subs": ["Dim1", "Dim2"]}, {},
-            _subscript_dict
+            {"expr": "1, 2; 3, 4; 5, 6;", "subs": ["Dim1", "Dim2"],
+             "merge_subs": ["Dim1", "Dim2"], "py_name": "var"},
+            {}, _subscript_dict,
         )
         string = element[0]["py_expr"]
         a = eval(string)
@@ -626,15 +660,129 @@ def test_subscript_3d_depth(self):
         self.assertEqual(a.loc[{"Dim1": "A", "Dim2": "D"}], 1)
         self.assertEqual(a.loc[{"Dim1": "B", "Dim2": "E"}], 4)
 
+    def test_subscript_builder(self):
+        """
+        Testing how subscripts are translated when we have common subscript
+        ranges.
+        """
+        from pysd.py_backend.vensim.vensim2py import parse_general_expression,\
+            parse_lookup_expression
+
+        _subscript_dict = {
+            "Dim1": ["A", "B", "C"], "Dim2": ["B", "C"], "Dim3": ["B", "C"]
+        }
+
+        # case 1: subscript of the expr is in the final range, which is a
+        # subrange of a greater range
+        element = parse_general_expression(
+            {"py_name": "var1", "subs": ["B"], "real_name": "var1", "eqn": "",
+             "expr": "GET DIRECT CONSTANTS('input.xlsx', 'Sheet1', 'C20')",
+             "merge_subs": ["Dim2"]},
+            {},
+            _subscript_dict
+        )
+        self.assertIn(
+            "'Dim2': ['B']", element[1][0]['py_expr'])
+
+        # case 1b: subscript of the expr is in the final range, which is a
+        # subrange of a greater range
+        element = parse_lookup_expression(
+            {"py_name": "var1b", "subs": ["B"],
+             "real_name": "var1b", "eqn": "",
+             "expr": "(GET DIRECT LOOKUPS('input.xlsx', 'Sheet1',"
+                     " '19', 'C20'))",
+             "merge_subs": ["Dim2"]},
+            _subscript_dict,
+        )
+        self.assertIn(
+            "'Dim2': ['B']", element[1][0]['py_expr'])
+
+        # case 2: subscript of the expr is a subscript subrange equal to the
+        # final range, which is a subrange of a greater range
+        element = parse_general_expression(
+            {"py_name": "var2", "subs": ["Dim2"],
+             "real_name": "var2", "eqn": "",
+             "expr": "GET DIRECT CONSTANTS('input.xlsx', 'Sheet1', 'C20')",
+             "merge_subs": ["Dim2"]},
+            {},
+            _subscript_dict
+        )
+        self.assertIn(
+            "'Dim2': _subscript_dict['Dim2']", element[1][0]['py_expr'])
+
+        # case 3: subscript of the expr is a subscript subrange equal to the
+        # final range, which is a subrange of a greater range, but there is
+        # a similar subrange before
+        element = parse_general_expression(
+            {"py_name": "var3", "subs": ["B"], "real_name": "var3", "eqn": "",
+             "expr": "GET DIRECT CONSTANTS('input.xlsx', 'Sheet1', 'C20')",
+             "merge_subs": ["Dim3"]},
+            {},
+            _subscript_dict
+        )
+        self.assertIn(
+            "'Dim3': ['B']", element[1][0]['py_expr'])
+
+        # case 4: subscript of the expr is a subscript subrange and the final
+        # subscript is a greater range
+        element = parse_general_expression(
+            {"py_name": "var4", "subs": ["Dim2"],
+             "real_name": "var4", "eqn": "",
+             "expr": "GET DIRECT CONSTANTS('input.xlsx', 'Sheet1', 'C20')",
+             "merge_subs": ["Dim1"]},
+            {},
+            _subscript_dict,
+        )
+        self.assertIn(
+            "'Dim1': _subscript_dict['Dim2']", element[1][0]['py_expr'])
+
+        # case 4b: subscript of the expr is a subscript subrange and the final
+        # subscript is a greater range
+        element = parse_general_expression(
+            {"py_name": "var4b", "subs": ["Dim2"],
+             "real_name": "var4b", "eqn": "",
+             "expr": "GET DIRECT DATA('input.xlsx', 'Sheet1', '19', 'C20')",
+             "keyword": None, "merge_subs": ["Dim1"]},
+            {},
+            _subscript_dict
+        )
+        self.assertIn(
+            "'Dim1': _subscript_dict['Dim2']", element[1][0]['py_expr'])
+
+        # case 4c: subscript of the expr is a subscript subrange and the final
+        # subscript is a greater range
+        element = parse_general_expression(
+            {"py_name": "var4c", "subs": ["Dim2"],
+             "real_name": "var4c", "eqn": "",
+             "expr": "GET DIRECT LOOKUPS('input.xlsx', 'Sheet1',"
+             " '19', 'C20')", "merge_subs": ["Dim1"]},
+            {},
+            _subscript_dict
+        )
+        self.assertIn(
+            "'Dim1': _subscript_dict['Dim2']", element[1][0]['py_expr'])
+
+        # case 4d: subscript of the expr is a subscript subrange and the final
+        # subscript is a greater range
+        element = parse_lookup_expression(
+            {"py_name": "var4d", "subs": ["Dim2"],
+             "real_name": "var4d", "eqn": "",
+             "expr": "(GET DIRECT LOOKUPS('input.xlsx', 'Sheet1',"
+                     " '19', 'C20'))", "merge_subs": ["Dim1"]},
+            _subscript_dict
+        )
+        self.assertIn(
+            "'Dim1': _subscript_dict['Dim2']", element[1][0]['py_expr'])
+
     def test_subscript_reference(self):
         from pysd.py_backend.vensim.vensim2py import parse_general_expression
 
         res = parse_general_expression(
-            {"expr": "Var A[Dim1, Dim2]"},
+            {"expr": "Var A[Dim1, Dim2]", "real_name": "Var2", "eqn": ""},
             {"Var A": "var_a"},
             {"Dim1": ["A", "B"], "Dim2": ["C", "D", "E"]},
             None,
-            {"var_a": ["Dim1", "Dim2"]},
+            {"var_a": ["Dim1", "Dim2"]}
         )
 
         self.assertEqual(res[0]["py_expr"], "var_a()")
@@ -712,7 +860,8 @@ def test_invert_matrix(self):
             {
                 "expr": "INVERT MATRIX(A, 3)",
                 "real_name": "A1",
-                "py_name": "a",
+                "py_name": "a1",
+                "merge_subs": ["dim1", "dim2"]
             },
             {
                 "A": "a",
@@ -720,10 +869,6 @@ def test_invert_matrix(self):
             },
             subscript_dict={
                 "dim1": ["a", "b", "c"], "dim2": ["a", "b", "c"]
-            },
-            elements_subs_dict={
-                "a1": ["dim1", "dim2"],
-                "a": ["dim1", "dim2"]
             }
         )
 
@@ -737,15 +882,13 @@ def test_subscript_logicals(self):
                 "expr": "IF THEN ELSE(dim1=dim2, 5, 0)",
                 "real_name": "A",
                 "py_name": "a",
+                "merge_subs": ["dim1", "dim2"]
             },
             {
                 "A": "a",
             },
             subscript_dict={
                 "dim1": ["a", "b", "c"], "dim2": ["a", "b", "c"]
-            },
-            elements_subs_dict={
-                "a": ["dim1", "dim2"]
             }
         )
 
@@ -766,13 +909,16 @@ def test_incomplete_expression(self):
                     "expr": "A FUNCTION OF(Unspecified Eqn,Var A,Var B)",
                     "real_name": "Incomplete Func",
                     "py_name": "incomplete_func",
+                    "eqn": "Incomplete Func = A FUNCTION OF(Unspecified "
+                           + "Eqn,Var A,Var B)",
+                    "subs": [],
+                    "merge_subs": []
                 },
                 {
                     "Unspecified Eqn": "unspecified_eqn",
                     "Var A": "var_a",
                     "Var B": "var_b",
-                },
-                elements_subs_dict={"incomplete_func": []},
+                }
             )
             self.assertEqual(len(w), 1)
             self.assertTrue(