diff --git a/src/main/python/systemds/operator/algorithm/builtin/pageRank.py b/src/main/python/systemds/operator/algorithm/builtin/pageRank.py
index 5e03e9dd93d..d1f037b9353 100644
--- a/src/main/python/systemds/operator/algorithm/builtin/pageRank.py
+++ b/src/main/python/systemds/operator/algorithm/builtin/pageRank.py
@@ -30,9 +30,6 @@
 
 
 def pageRank(G: Matrix,
-             p: Matrix,
-             e: Matrix,
-             u: Matrix,
              **kwargs: Dict[str, VALID_INPUT_TYPES]):
     """
      DML builtin method for PageRank algorithm (power iterations)
@@ -41,14 +38,16 @@ def pageRank(G: Matrix,
     
     :param G: Input Matrix
     :param p: initial page rank vector (number of nodes), e.g., rand intialized
+        default rand initialized with seed
     :param e: additional customization, default vector of ones
-    :param u: personalization vector (number of nodes)
+    :param u: personalization vector (number of nodes), default vector of ones
     :param alpha: teleport probability
     :param max_iter: maximum number of iterations
+    :param seed: seed for default rand initialization of page rank vector
     :return: computed pagerank
     """
 
-    params_dict = {'G': G, 'p': p, 'e': e, 'u': u}
+    params_dict = {'G': G}
     params_dict.update(kwargs)
     return Matrix(G.sds_context,
         'pageRank',
diff --git a/src/main/python/systemds/operator/nn/__init__.py b/src/main/python/systemds/operator/nn/__init__.py
new file mode 100644
index 00000000000..e66abb4646f
--- /dev/null
+++ b/src/main/python/systemds/operator/nn/__init__.py
@@ -0,0 +1,20 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
diff --git a/src/main/python/systemds/operator/nn/affine.py b/src/main/python/systemds/operator/nn/affine.py
new file mode 100644
index 00000000000..44c67d1eda6
--- /dev/null
+++ b/src/main/python/systemds/operator/nn/affine.py
@@ -0,0 +1,114 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
+import os
+
+from systemds.context import SystemDSContext
+from systemds.operator import Matrix, Source, MultiReturn
+from systemds.utils.helpers import get_path_to_script_layers
+
+
+class Affine:
+    _source: Source = None
+    weight: Matrix
+    bias: Matrix
+
+    def __new__(cls, *args, **kwargs):
+        return super().__new__(cls)
+
+    def __init__(self, sds_context: SystemDSContext, d, m, seed=-1):
+        """
+        sds_context: The systemdsContext to construct the layer inside of
+        d: The number of features that are input to the affine layer
+        m: The number of neurons that are contained in the layer, 
+            and the number of features output
+        """
+        Affine._create_source(sds_context)
+
+        # bypassing overload limitation in python
+        self.forward = self._instance_forward
+        self.backward = self._instance_backward
+
+        # init weight and bias
+        self.weight = Matrix(sds_context, '')
+        self.bias = Matrix(sds_context, '')
+        params_dict = {'D': d, 'M': m, 'seed': seed}
+        out = [self.weight, self.bias]
+        op = MultiReturn(sds_context, "affine::init", output_nodes=out, named_input_nodes=params_dict)
+        self.weight._unnamed_input_nodes = [op]
+        self.bias._unnamed_input_nodes = [op]
+        op._source_node = self._source
+
+    @staticmethod
+    def forward(X: Matrix, W: Matrix, b: Matrix):
+        """
+        X: An input matrix
+        W: The hidden weights for the affine layer
+        b: The bias added in the output.
+        return out: An output matrix.
+        """
+        Affine._create_source(X.sds_context)
+        return Affine._source.forward(X, W, b)
+
+    @staticmethod
+    def backward(dout:Matrix, X: Matrix, W: Matrix, b: Matrix):
+        """
+        dout: The gradient of the output, passed from the upstream
+        X: The input matrix of this layer
+        W: The hidden weights for the affine layer
+        b: The bias added in the output
+        return dX, dW, db: The gradients of: input X, weights and bias.
+        """
+        sds = X.sds_context
+        Affine._create_source(sds)
+        params_dict = {'dout': dout, 'X': X, 'W': W, 'b': b}
+        dX = Matrix(sds, '')
+        dW = Matrix(sds, '')
+        db = Matrix(sds, '')
+        out = [dX, dW, db]
+        op = MultiReturn(sds, "affine::backward", output_nodes=out, named_input_nodes=params_dict)
+        dX._unnamed_input_nodes = [op]
+        dW._unnamed_input_nodes = [op]
+        db._unnamed_input_nodes = [op]
+        op._source_node = Affine._source
+        return op
+
+    def _instance_forward(self, X: Matrix):
+        """
+        X: The input matrix
+        return out: The output matrix
+        """
+        self._X = X
+        return Affine.forward(X, self.weight, self.bias)
+
+    def _instance_backward(self, dout: Matrix, X: Matrix):
+        """
+        dout: The gradient of the output, passed from the upstream layer
+        X: The input to this layer.
+        return dX, dW,db: gradient of input, weights and bias, respectively
+        """
+        return Affine.backward(dout, X, self.weight, self.bias)
+
+    @staticmethod
+    def _create_source(sds: SystemDSContext):
+        if Affine._source is None or Affine._source.sds_context != sds:
+            path = get_path_to_script_layers()
+            path = os.path.join(path, "affine.dml")
+            Affine._source = sds.source(path, "affine")
diff --git a/src/main/python/systemds/operator/nn/relu.py b/src/main/python/systemds/operator/nn/relu.py
new file mode 100644
index 00000000000..99833e6d86d
--- /dev/null
+++ b/src/main/python/systemds/operator/nn/relu.py
@@ -0,0 +1,68 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
+import os.path
+
+from systemds.context import SystemDSContext
+from systemds.operator import Matrix, Source
+from systemds.utils.helpers import get_path_to_script_layers
+
+
+class ReLU:
+    _source: Source = None
+
+    def __init__(self, sds: SystemDSContext):
+        ReLU._create_source(sds)
+        self.forward = self._instance_forward
+        self.backward = self._instance_backward
+
+    @staticmethod
+    def forward(X: Matrix):
+        """
+        X: input matrix
+        return out: output matrix
+        """
+        ReLU._create_source(X.sds_context)
+        return ReLU._source.forward(X)
+
+    @staticmethod
+    def backward(dout: Matrix, X: Matrix):
+        """
+        dout: gradient of output, passed from the upstream
+        X: input matrix
+        return dX: gradient of input
+        """
+        ReLU._create_source(dout.sds_context)
+        return ReLU._source.backward(dout, X)
+
+    def _instance_forward(self, X: Matrix):
+        self._X = X
+        return ReLU.forward(X)
+
+    def _instance_backward(self, dout: Matrix, X: Matrix):
+        return ReLU.backward(dout, X)
+
+    @staticmethod
+    def _create_source(sds: SystemDSContext):
+        if ReLU._source is None or ReLU._source.sds_context != sds:
+            path = get_path_to_script_layers()
+            path = os.path.join(path, "relu.dml")
+            ReLU._source = sds.source(path, "relu")
+
diff --git a/src/main/python/systemds/operator/nodes/source.py b/src/main/python/systemds/operator/nodes/source.py
index 7191fc8208e..388ae9fd376 100644
--- a/src/main/python/systemds/operator/nodes/source.py
+++ b/src/main/python/systemds/operator/nodes/source.py
@@ -21,6 +21,8 @@
 
 __all__ = ["Source"]
 
+import platform
+
 from types import MethodType
 from typing import TYPE_CHECKING, Dict, Iterable, Sequence
 
@@ -142,6 +144,7 @@ def __init__(self, sds_context, path: str, name: str):
             func = f.get_func(sds_context, name)
             setattr(self, f._name, MethodType(func, self))
 
+
     def __parse_functions_from_script(self, path: str) -> Iterable[Func]:
         lines = self.__parse_lines_with_filter(path)
         functions = []
@@ -162,10 +165,18 @@ def __parse_lines_with_filter(self, path: str) -> Iterable[str]:
         lines = []
         with open(path) as file:
             insideBracket = 0
+            insideComment = False
             for l in file.readlines():
                 ls = l.strip()
                 if len(ls) == 0 or ls[0] == '#':
                     continue
+                elif insideComment:
+                    if ls.endswith("*/"):
+                        insideComment = False
+                        continue
+                elif ls.startswith("/*"):
+                    insideComment = True
+                    continue
                 elif insideBracket > 0:
                     for c in ls:
                         if c == '{':
@@ -193,7 +204,11 @@ def __parse_lines_with_filter(self, path: str) -> Iterable[str]:
         return filtered_lines
 
     def code_line(self, var_name: str, unnamed_input_vars: Sequence[str], named_input_vars: Dict[str, str]) -> str:
-        line = f'source({self.operation}) as { self.__name}'
+        if platform.system() == 'Windows':
+            source_path = self.operation.replace("\\","\\\\")
+        else:
+            source_path = self.operation
+        line = f'source({source_path}) as { self.__name}'
         return line
 
     def compute(self, verbose: bool = False, lineage: bool = False):
diff --git a/src/main/python/systemds/utils/helpers.py b/src/main/python/systemds/utils/helpers.py
index b25ac65735c..4f55a006233 100644
--- a/src/main/python/systemds/utils/helpers.py
+++ b/src/main/python/systemds/utils/helpers.py
@@ -64,10 +64,10 @@ def get_slice_string(i):
             raise NotImplementedError("Not Implemented slice with dynamic end")
         else:
             # + 1 since R and systemDS is 1 indexed.
-            return f'{i.start+1}:{i.stop}'
+            return f'{i.start + 1}:{i.stop}'
     else:
         # + 1 since R and systemDS is 1 indexed.
-        sliceIns = i+1
+        sliceIns = i + 1
     return sliceIns
 
 
@@ -77,5 +77,12 @@ def check_is_empty_slice(i):
 
 def check_no_less_than_zero(i: list):
     for x in i:
-        if(x < 0):
+        if (x < 0):
             raise ValueError("Negative index not supported in systemds")
+
+
+def get_path_to_script_layers() -> str:
+    root = os.environ.get("SYSTEMDS_ROOT")
+    if root is None:
+        root = get_module_dir()
+    return os.path.join(root, "scripts", "nn", "layers")
diff --git a/src/main/python/tests/nn/__init__.py b/src/main/python/tests/nn/__init__.py
new file mode 100644
index 00000000000..e66abb4646f
--- /dev/null
+++ b/src/main/python/tests/nn/__init__.py
@@ -0,0 +1,20 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
diff --git a/src/main/python/tests/nn/neural_network.py b/src/main/python/tests/nn/neural_network.py
new file mode 100644
index 00000000000..e62236e9bcd
--- /dev/null
+++ b/src/main/python/tests/nn/neural_network.py
@@ -0,0 +1,89 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
+from systemds.context import SystemDSContext
+
+from systemds.operator.nn.affine import Affine
+from systemds.operator.nn.relu import ReLU
+from systemds.operator import Matrix, Source
+
+
+class NeuralNetwork:
+    _source: Source = None
+    _X: Matrix
+
+    def __init__(self, sds: SystemDSContext, dim: int):
+
+        # first hidden layer
+        self.affine1 = Affine(sds, dim, 128, seed=42)
+        self.w1, self.b1 = self.affine1.weight, self.affine1.bias
+        self.relu1 = ReLU(sds)
+
+        # second hidden layer
+        self.affine2 = Affine(sds, 128, 64, seed=42)
+        self.w2, self.b2 = self.affine2.weight, self.affine2.bias
+        self.relu2 = ReLU(sds)
+
+        # third hidden layer
+        self.affine3 = Affine(sds, 64, 32, seed=42)
+        self.relu3 = ReLU(sds)
+        self.w3, self.b3 = self.affine3.weight, self.affine3.bias
+
+        # output layer
+        self.affine4 = Affine(sds, 32, 2, seed=42)
+        self.w4, self.b4 = self.affine4.weight, self.affine4.bias
+
+    def forward_static_pass(self, X: Matrix) -> Matrix:
+        """
+        Compute forward pass through the network using static affine and relu calls
+        :param X: Input matrix
+        :return: Output matrix
+        """
+        X = self.affine1.forward(X)
+        X = self.relu1.forward(X)
+
+        X = self.affine2.forward(X)
+        X = self.relu2.forward(X)
+
+        X = self.affine3.forward(X)
+        X = self.relu3.forward(X)
+
+        X = self.affine4.forward(X)
+
+        return X
+
+    def forward_dynamic_pass(self, X: Matrix) -> Matrix:
+        """
+        Compute forward pass through the network using dynamic affine and relu calls
+        :param X: Input matrix
+        :return: Output matrix
+        """
+        X = Affine.forward(X, self.w1, self.b1)
+        X = ReLU.forward(X)
+
+        X = Affine.forward(X, self.w2, self.b2)
+        X = ReLU.forward(X)
+
+        X = Affine.forward(X, self.w3, self.b3)
+        X = ReLU.forward(X)
+
+        X = Affine.forward(X, self.w4, self.b4)
+
+        return X
diff --git a/src/main/python/tests/nn/test_affine.py b/src/main/python/tests/nn/test_affine.py
new file mode 100644
index 00000000000..955945b29c1
--- /dev/null
+++ b/src/main/python/tests/nn/test_affine.py
@@ -0,0 +1,163 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
+
+import unittest
+
+import numpy as np
+from numpy.testing import assert_almost_equal
+
+from systemds.context import SystemDSContext
+from systemds.script_building.script import DMLScript
+from systemds.operator.nn.affine import Affine
+
+dim = 6
+n = 5
+m = 6
+X = np.array([[9., 2., 5., 5., 9., 6.],
+              [0., 8., 8., 0., 5., 7.],
+              [2., 2., 6., 3., 4., 3.],
+              [3., 5., 2., 6., 6., 0.],
+              [3., 8., 5., 2., 5., 2.]])
+
+W = np.array([[8., 3., 7., 2., 0., 1.],
+              [6., 5., 1., 2., 6., 1.],
+              [2., 4., 7., 7., 6., 4.],
+              [3., 8., 9., 3., 5., 6.],
+              [3., 8., 0., 5., 7., 9.],
+              [7., 9., 7., 4., 5., 7.]])
+dout = np.array([[9., 5., 4., 0., 4., 1.],
+                 [1., 2., 2., 3., 3., 9.],
+                 [7., 4., 0., 8., 7., 0.],
+                 [8., 7., 0., 6., 0., 9.],
+                 [1., 6., 5., 8., 8., 9.]])
+
+
+class TestAffine(unittest.TestCase):
+    sds: SystemDSContext = None
+
+    @classmethod
+    def setUpClass(cls):
+        cls.sds = SystemDSContext()
+
+    @classmethod
+    def tearDownClass(cls):
+        cls.sds.close()
+
+    def test_init(self):
+        affine = Affine(self.sds, dim, m, 10)
+        w = affine.weight.compute()
+        self.assertEqual(len(w), 6)
+        self.assertEqual(len(w[0]), 6)
+
+    def test_forward(self):
+        Xm = self.sds.from_numpy(X)
+        Wm = self.sds.from_numpy(W)
+        bm = self.sds.full((1, 6), 0)
+
+        # test class method
+        affine = Affine(self.sds, dim, m, 10)
+        out = affine.forward(Xm).compute()
+        self.assertEqual(len(out), 5)
+        self.assertEqual(len(out[0]), 6)
+
+        # test static method
+        out = Affine.forward(Xm, Wm, bm).compute()
+        expected = np.matmul(X, W)
+        assert_almost_equal(out, expected)
+
+    def test_backward(self):
+        Xm = self.sds.from_numpy(X)
+        Wm = self.sds.from_numpy(W)
+        bm = self.sds.full((1, 6), 0)
+        doutm = self.sds.from_numpy(dout)
+
+        # test class method
+        affine = Affine(self.sds, dim, m, 10)
+        [dx, dw, db] = affine.backward(doutm, Xm).compute()
+        assert len(dx) == 5 and len(dx[0]) == 6
+        assert len(dw) == 6 and len(dx[0]) == 6
+        assert len(db) == 1 and len(db[0]) == 6
+
+        # test static method
+        res = Affine.backward(doutm, Xm, Wm, bm).compute()
+        assert len(res) == 3
+
+    def test_multiple_sourcing(self):
+        sds = SystemDSContext()
+        a1 = Affine(sds, dim, m, 10)
+        a2 = Affine(sds, m, 11, 10)
+
+        Xm = sds.from_numpy(X)
+        X1 = a1.forward(Xm)
+        X2 = a2.forward(X1)
+
+        scripts = DMLScript(sds)
+        scripts.build_code(X2)
+
+        self.assertEqual(1, self.count_sourcing(scripts.dml_script, layer_name="affine"))
+        sds.close()
+
+    def test_multiple_context(self):
+        # This test evaluate if multiple conflicting contexts work.
+        # It is not the 'optimal' nor the intended use 
+        # If it fails in the future, feel free to delete it.
+        
+        # two context
+        sds1 = SystemDSContext()
+        sds2 = SystemDSContext()
+        a1 = Affine(sds1, dim, m, 10)
+        a2 = Affine(sds2, m, 11, 10)
+
+        Xm = sds1.from_numpy(X)
+        X1 = a1.forward(Xm)
+        out_actual = a2.forward(X1).compute()
+
+        # one context
+        Xm = self.sds.from_numpy(X)
+        a1 = Affine(self.sds, dim, m, 10)
+        a2 = Affine(self.sds, m, 11, 10)
+
+        X1 = a1.forward(Xm)
+        out_expected = a2.forward(X1).compute()
+
+        assert_almost_equal(out_actual, out_expected)
+
+        sds1.close()
+        sds2.close()
+
+    def count_sourcing(self, script: str, layer_name: str):
+        """
+        Count the number of times the dml script is being sourced
+        i.e. count the number of occurrences of lines like
+        'source(...) as affine' in the dml script
+
+        :param script: the sourced dml script text
+        :param layer_name: example: "affine", "relu"
+        :return:
+        """
+        return len([
+            line for line in script.split("\n")
+            if all([line.startswith("source"), line.endswith(layer_name)])
+        ])
+
+
+if __name__ == '__main__':
+    unittest.main()
diff --git a/src/main/python/tests/nn/test_neural_network.py b/src/main/python/tests/nn/test_neural_network.py
new file mode 100644
index 00000000000..1cf3c5dd332
--- /dev/null
+++ b/src/main/python/tests/nn/test_neural_network.py
@@ -0,0 +1,94 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
+import unittest
+import numpy as np
+
+from systemds.context import SystemDSContext
+from tests.nn.neural_network import NeuralNetwork
+from systemds.script_building.script import DMLScript
+
+# Seed for the input matrix
+np.random.seed(42)
+
+
+class TestNeuralNetwork(unittest.TestCase):
+    sds: SystemDSContext = None
+
+    @classmethod
+    def setUpClass(cls):
+        cls.sds = SystemDSContext()
+        cls.X = np.random.rand(6, 1)
+        cls.exp_out = np.array([
+            -0.37768756, -0.47785831, -0.95870362,
+            -1.21297214, -0.73814523, -0.933917,
+            -0.60368929, -0.76380049, -0.15732974,
+            -0.19905692, -0.15730542, -0.19902615
+        ])
+
+    @classmethod
+    def tearDownClass(cls):
+        cls.sds.close()
+
+    def test_forward_pass(self):
+
+        Xm = self.sds.from_numpy(self.X)
+        nn = NeuralNetwork(self.sds, dim=1)
+        # test forward pass through the network using static calls
+        static_out = nn.forward_static_pass(Xm).compute().flatten()
+
+        self.assertTrue(np.allclose(static_out, self.exp_out))
+
+        # test forward pass through the network using dynamic calls
+        dynamic_out = nn.forward_dynamic_pass(Xm).compute().flatten()
+        self.assertTrue(np.allclose(dynamic_out,self.exp_out))
+
+    def test_multiple_sourcing(self):
+        sds = SystemDSContext()
+        Xm = sds.from_numpy(self.X)
+        nn = NeuralNetwork(sds, dim=1)
+
+        # test for verifying that affine and relu are each being sourced exactly once
+        network_out = nn.forward_static_pass(Xm)
+        scripts = DMLScript(sds)
+        scripts.build_code(network_out)
+
+        self.assertEqual(1, self.count_sourcing(scripts.dml_script, layer_name="affine"))
+        self.assertEqual(1, self.count_sourcing(scripts.dml_script, layer_name="relu"))
+        sds.close()
+
+    def count_sourcing(self, script: str, layer_name: str):
+        """
+        Count the number of times the dml script is being sourced
+        i.e. count the number of occurrences of lines like
+        'source(...) as relu' in the dml script
+
+        :param script: the sourced dml script text
+        :param layer_name: example: "affine", "relu"
+        :return:
+        """
+        return len([
+            line for line in script.split("\n")
+            if all([line.startswith("source"), line.endswith(layer_name)])
+        ])
+
+
+if __name__ == '__main__':
+    unittest.main()
diff --git a/src/main/python/tests/nn/test_relu.py b/src/main/python/tests/nn/test_relu.py
new file mode 100644
index 00000000000..06839ce494b
--- /dev/null
+++ b/src/main/python/tests/nn/test_relu.py
@@ -0,0 +1,105 @@
+# -------------------------------------------------------------
+#
+# Licensed to the Apache Software Foundation (ASF) under one
+# or more contributor license agreements.  See the NOTICE file
+# distributed with this work for additional information
+# regarding copyright ownership.  The ASF licenses this file
+# to you under the Apache License, Version 2.0 (the
+# "License"); you may not use this file except in compliance
+# with the License.  You may obtain a copy of the License at
+#
+#   http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing,
+# software distributed under the License is distributed on an
+# "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
+# KIND, either express or implied.  See the License for the
+# specific language governing permissions and limitations
+# under the License.
+#
+# -------------------------------------------------------------
+
+import unittest
+
+import numpy as np
+
+from systemds.context import SystemDSContext
+from systemds.script_building.script import DMLScript
+from systemds.operator.nn.relu import ReLU
+
+
+class TestRelu(unittest.TestCase):
+    sds: SystemDSContext = None
+
+    @classmethod
+    def setUpClass(cls):
+        cls.sds = SystemDSContext()
+
+    @classmethod
+    def tearDownClass(cls):
+        cls.sds.close()
+
+    def test_forward(self):
+        X = np.array([0, -1, -2, 2, 3, -5])
+        relu = ReLU(self.sds)
+        # forward
+        Xm = self.sds.from_numpy(X)
+        out = relu.forward(Xm).compute().flatten()
+        expected = np.array([0, 0, 0, 2, 3, 0])
+        self.assertTrue(np.allclose(out, expected))
+
+        # test static
+        sout = ReLU.forward(Xm).compute().flatten()
+        self.assertTrue(np.allclose(sout, expected))
+
+    def test_backward(self):
+        X = np.array([0, -1, -2, 2, 3, -5])
+        dout = np.array([0, 1, 2, 3, 4, 5])
+        relu = ReLU(self.sds)
+        # forward
+        Xm = self.sds.from_numpy(X)
+        out = relu.forward(Xm)
+        # backward
+        doutm = self.sds.from_numpy(dout)
+        dx = relu.backward(doutm, Xm).compute().flatten()
+        expected = np.array([0, 0, 0, 3, 4, 0], dtype=np.double)
+        self.assertTrue(np.allclose(dx, expected))
+
+        # test static
+        sdx = ReLU.backward(doutm, Xm).compute().flatten()
+        self.assertTrue(np.allclose(sdx, expected))
+
+    def test_multiple_sourcing(self):
+        sds = SystemDSContext()
+        X = np.array([0, -1, -2, 2, 3, -5])
+        r1 = ReLU(sds)
+        r2 = ReLU(sds)
+
+        Xm = sds.from_numpy(X)
+        X1 = r1.forward(Xm)
+        X2 = r2.forward(X1)
+
+        scripts = DMLScript(sds)
+        scripts.build_code(X2)
+
+        self.assertEqual(1,self.count_sourcing(scripts.dml_script, layer_name="relu"))
+        sds.close()
+
+    def count_sourcing(self, script: str, layer_name: str):
+        """
+        Count the number of times the dml script is being sourced
+        i.e. count the number of occurrences of lines like
+        'source(...) as relu' in the dml script
+
+        :param script: the sourced dml script text
+        :param layer_name: example: "affine", "relu"
+        :return:
+        """
+        return len([
+            line for line in script.split("\n")
+            if all([line.startswith("source"), line.endswith(layer_name)])
+        ])
+
+
+if __name__ == '__main__':
+    unittest.main()