diff --git a/.github/environment.yml b/.github/environment.yml
index 89392b5..5543a65 100644
--- a/.github/environment.yml
+++ b/.github/environment.yml
@@ -4,7 +4,7 @@ channels:
   - conda-forge
 
 dependencies:
-  - python==3.8
+  - python=PYTHON_VERSION
   - pip
   - oommf
   - pip:
diff --git a/.github/workflows/workflow-gpu.yml b/.github/workflows/workflow-gpu.yml
index 0a1926a..2d4cba9 100644
--- a/.github/workflows/workflow-gpu.yml
+++ b/.github/workflows/workflow-gpu.yml
@@ -64,6 +64,9 @@ jobs:
     defaults:
       run:
         shell: bash -i {0}
+    strategy:
+      matrix:
+        python-version: ["3.8", "3.10"]
 
     steps:
     - name: Initialisation
@@ -73,6 +76,10 @@ jobs:
         pwd
         echo $PATH
         which conda
+
+    - name: Prepare env file
+      run: sed -i "s/PYTHON_VERSION/${{ matrix.python-version }}/g" .github/environment.yml
+
     - name: Create conda env
       run: conda env create -f ".github/environment.yml"
 
diff --git a/README.md b/README.md
index 03742c0..101d8bd 100644
--- a/README.md
+++ b/README.md
@@ -83,3 +83,4 @@ Licensed under the BSD 3-Clause "New" or "Revised" License. For details, please
 - [OpenDreamKit](http://opendreamkit.org/) – Horizon 2020 European Research Infrastructure project (676541)
 
 - EPSRC Programme Grant on [Skyrmionics](http://www.skyrmions.ac.uk) (EP/N032128/1)
+- 
diff --git a/mumax3c/scripts/driver.py b/mumax3c/scripts/driver.py
index a8f67d5..a4c1ccc 100644
--- a/mumax3c/scripts/driver.py
+++ b/mumax3c/scripts/driver.py
@@ -11,6 +11,7 @@ def driver_script(driver, system, compute=None, ovf_format="bin4", **kwargs):
     mx3 = "tableadd(E_total)\n"
     mx3 += "tableadd(dt)\n"
     mx3 += "tableadd(maxtorque)\n"
+
     if isinstance(driver, mc.MinDriver):
         for attr, value in driver:
             if attr != "evolver":
@@ -34,6 +35,9 @@ def driver_script(driver, system, compute=None, ovf_format="bin4", **kwargs):
         mx3 += "tablesave()\n\n"
 
     if isinstance(driver, mc.TimeDriver):
+        # Need temperature only here
+        if system.T > 0:
+            mx3 += f"Temp = {system.T}\n"
         # Extract dynamics equation parameters.
         gamma0 = (
             precession[0].gamma0
diff --git a/mumax3c/tests/test_micromagnetictests.py b/mumax3c/tests/test_micromagnetictests.py
index 41a9914..9b34d11 100644
--- a/mumax3c/tests/test_micromagnetictests.py
+++ b/mumax3c/tests/test_micromagnetictests.py
@@ -1 +1,32 @@
+import discretisedfield as df
+import micromagneticmodel as mm
+import numpy as np
 from micromagnetictests.calculatortests import *  # noqa: F401,F403
+
+import mumax3c as mc
+
+
+def test_temperature():
+    p1 = (0, 0, 0)
+    p2 = (5e-9, 5e-9, 5e-9)
+    n = (2, 2, 2)
+    Ms = 1e6
+    A = 1e-13
+    region = df.Region(p1=p1, p2=p2)
+    mesh = df.Mesh(region=region, n=n)
+
+    system = mm.System(name="mumax3_temperature")
+    system.energy = mm.Exchange(A=A)
+    system.dynamics = mm.Precession(gamma0=mm.consts.gamma0) + mm.Damping(alpha=1)
+    system.m = df.Field(mesh, nvdim=3, value=(0, 0.1, 1), norm=Ms)
+    system.T = 1000  # Set a high temperature so it is almost random
+    assert f"Temp = {system.T}" in mc.scripts.driver_script(
+        mc.TimeDriver(), system, compute=None, t=1, n=1
+    )
+
+    td = mc.TimeDriver()
+    td.drive(system, t=0.2e-9, n=1)
+
+    # If random then there should be no overall direction
+    assert np.linalg.norm(system.m.orientation.mean()) < 0.5
+    mc.delete(system)