diff --git a/models/bonemodel2.py b/models/bonemodel2.py
index 187fb64ee6adb3e3e6ca18a0d315a144a36c99c3..5a08f4d2c8b2849565a24b3c144f4dbbc8a0a31c 100644
--- a/models/bonemodel2.py
+++ b/models/bonemodel2.py
@@ -47,28 +47,11 @@ def parms(d={}):
         'axis_pt1': ['x', 'y', 'z'],
         'axis_pt2': ['x', 'y', 'z']
     }
-    # rigid cylinder (used if p['food']=='cylinder')
-    # R = 20
-    # p['cylinder'] = {
-    #     'origin': [-R-8, 0, -150],
-    #     'radius': R,
-    #     'width': 100,
-    #     'friction': 0.1,
-    #     'penalty': 1e3
-    # }
 
     # material properties
-    # p['density'] = 1.850e-9  # [T/mm³] - bone: 1.850 kg/l
     p['Young'] = 17000.      # [MPa] elastic modulus - bone: 17-20 GPa
     p['Poisson'] = 0.3       # [-] Poisson's ratio
 
-    # numerical parameters
-    # p['tolNR'] = 1e-6        # [-] equilibrium tolerance
-    # p['dt0'] = 1.0           # [s] time step size
-
-    # gmsh toolbox
-    # p['use_gmshOld'] = False  # use old gmsh interface
-
     p.update(d)
     return p
 
@@ -130,6 +113,12 @@ def solve(p={}):
         # store data in a structure
         mgroups[name] = MuscleGroup(nodes, tris, ntags, loads)
 
+        # save loads to file
+        with open(f'loads_{name}.tsv','w') as f:
+            f.write('node#\tfx\tfy\tfz\n')
+            for i,load in enumerate(loads):
+                f.write(f"{int(ntags[i])}\t{load.x[0]}\t{load.x[1]}\t{load.x[2]}\n")   
+
     if p['food'] == 'fixteeth':
         # teeth are fixed along chosen directions
         create_group(p['contact_pts'], nods_no,