From 7a9243edffc8af7ad0122d5b23bcfcab923233e6 Mon Sep 17 00:00:00 2001
From: Romain Boman <r.boman@uliege.be>
Date: Mon, 29 Mar 2021 11:53:08 +0200
Subject: [PATCH] add fixations

---
 models/dolicorhynchops_10k.py         | 13 ++++++++++++-
 models/dolicorhynchops_10k_contact.py | 14 +++++++++++++-
 models/dolicorhynchops_149k.py        | 11 +++++++++++
 models/mandiblemodel.py               | 16 +++++++++++-----
 4 files changed, 47 insertions(+), 7 deletions(-)

diff --git a/models/dolicorhynchops_10k.py b/models/dolicorhynchops_10k.py
index f8ce731..fb8dcbf 100644
--- a/models/dolicorhynchops_10k.py
+++ b/models/dolicorhynchops_10k.py
@@ -16,7 +16,7 @@ def parms(d={}):
             'file': f'{path}/Lmuscle.ply',
             'force': 100.,
             'focalpt': f'{path}/LmuscleF.off',
-            'method': 'T'                              # 'U', 'T', 'T+N'
+            'method': 'U'                              # 'U', 'T', 'T+N'
         },
         {
             'file': f'{path}/Rmuscle.off',
@@ -25,6 +25,17 @@ def parms(d={}):
             'method': 'T'
         }
     ]
+    p['fixations'] = {  
+        'teeth': ['x'],
+        'LTMJ': ['x','y','z'],
+        'RTMJ': ['x','z']
+    }
+
+    # material properties
+    p['density'] = 1.850e-9  # [T/mm³] - bone: 1.850 kg/l
+    p['Young'] = 20000.      # [MPa] elastic modulus - bone: 17-20 GPa
+    p['Poisson'] = 0.3       # [-] Poisson's ratio
+
     p.update(d)
     return p
 
diff --git a/models/dolicorhynchops_10k_contact.py b/models/dolicorhynchops_10k_contact.py
index 2aa528c..0ed0880 100644
--- a/models/dolicorhynchops_10k_contact.py
+++ b/models/dolicorhynchops_10k_contact.py
@@ -9,7 +9,7 @@ def parms(d={}):
     path = 'dolicorhynchops/10k'
     p['mandible'] = f'{path}/mandible.ply'
     # p['teeth'] = f'{path}/teeth.off'
-    p['LTMJ'] = f'{path}/LTMJ.off'
+    p['LTMJ'] = [0.1458893, -73.13895, 227.3909],
     p['RTMJ'] = f'{path}/RTMJ.off'
     p['muscles'] = [
         {
@@ -25,6 +25,11 @@ def parms(d={}):
             'method': 'T'
         }
     ]
+    p['fixations'] = {  
+        'teeth': ['x','y','z'],
+        'LTMJ': ['x','y','z'],
+        'RTMJ': ['x','y','z']
+    }
     p['food'] = 'cylinder'
     R = 10
     p['cylinder'] = {
@@ -37,6 +42,13 @@ def parms(d={}):
     p['tolNR'] = 1e-4        # [-] equilibrium tolerance
     p['dt0'] = 2e-2          # [s] time step size
 
+
+    # material properties
+    p['density'] = 1.850e-9  # [T/mm³] - bone: 1.850 kg/l
+    p['Young'] = 20000.      # [MPa] elastic modulus - bone: 17-20 GPa
+    p['Poisson'] = 0.3       # [-] Poisson's ratio
+
+
     p.update(d)
     return p
 
diff --git a/models/dolicorhynchops_149k.py b/models/dolicorhynchops_149k.py
index c47661f..53c839e 100644
--- a/models/dolicorhynchops_149k.py
+++ b/models/dolicorhynchops_149k.py
@@ -25,6 +25,17 @@ def parms(d={}):
             'method': 'T+N'
         }
     ]
+    p['fixations'] = {  
+        'teeth': ['x'],
+        'LTMJ': ['x','y','z'],
+        'RTMJ': ['x','z']
+    }
+
+    # material properties
+    p['density'] = 1.850e-9  # [T/mm³] - bone: 1.850 kg/l
+    p['Young'] = 20000.      # [MPa] elastic modulus - bone: 17-20 GPa
+    p['Poisson'] = 0.3       # [-] Poisson's ratio
+
     p.update(d)
     return p
 
diff --git a/models/mandiblemodel.py b/models/mandiblemodel.py
index 9b84dd6..e1b7401 100644
--- a/models/mandiblemodel.py
+++ b/models/mandiblemodel.py
@@ -26,7 +26,11 @@ def parms(d={}):
         { 'file': f'{path}/Rmuscle.off', 'force': 100., 'focalpt':f'{path}/RmuscleF.off', 'method':'T'}
     ]
     p['food'] = 'fixteeth'  # type of food: 'fixteeth' or 'cylinder'
-
+    p['fixations'] = {  
+        'teeth': ['x','y','z'],  # constrained degrees of freedom
+        'LTMJ': ['x','y','z'],
+        'RTMJ': ['x','y','z']
+    }
     # rigid cylinder (used if p['food']=='cylinder') 
     R = 20
     p['cylinder'] = {
@@ -123,10 +127,12 @@ def getMetafor(p={}):
     fct.setData(0.0, 0.0)
     fct.setData(1.0, 1.0)
 
+    txt2key = { 'x': TX, 'y':TY, 'z':TZ }
+
     # axis of rotation
     for gname in ['LTMJ', 'RTMJ']:
-        for d in [ TX, TY, TZ ]:
-            domain.getLoadingSet().define(groups[gname], Field1D(d,RE), 0.0)
+        for d in p['fixations'][gname]:
+            domain.getLoadingSet().define(groups[gname], Field1D(txt2key[d],RE), 0.0)
 
     mshpts = geometry.getMesh().getPointSet()
     for name, gr in mgroups.items():
@@ -162,8 +168,8 @@ def getMetafor(p={}):
     elif p['food']=='fixteeth':
         # teeth are fixed along X
         create_group(p['teeth'], nods_no, nods_pos, domain, groups, 'teeth')
-        for d in [ TX ]:
-            domain.getLoadingSet().define(groups['teeth'], Field1D(d,RE), 0.0)
+        for d in p['fixations']['teeth']:
+            domain.getLoadingSet().define(groups['teeth'], Field1D(txt2key[d],RE), 0.0)
     else:
         # p['food'] is a filename
         ply2mtf(domain, p['food'])
-- 
GitLab