From 6cfd931bc450aa20326c16fdc0d528077b125e79 Mon Sep 17 00:00:00 2001 From: Romain Boman <romain.boman@gmail.com> Date: Mon, 6 May 2024 15:00:29 +0200 Subject: [PATCH] add new way of adding fixations --- models/bonemodel2.py | 71 +++++++++---------- models/dolicorhynchops/dolicorhynchops_10k.py | 27 ++++--- 2 files changed, 52 insertions(+), 46 deletions(-) diff --git a/models/bonemodel2.py b/models/bonemodel2.py index bd7cd8f..3f2634b 100644 --- a/models/bonemodel2.py +++ b/models/bonemodel2.py @@ -18,16 +18,11 @@ def parms(d={}): """ p = {} - # Geomagic/MeshLab input files + # Geomagic/MeshLab input folder path = 'dolicorhynchops/10k' + # mandible mesh (.ply/.geo/.msh) p['bone'] = f'{path}/mandible.ply' - # one or several vertices (.ply/.off) (used if p['food']=='fixteeth') - p['contact_pts'] = f'{path}/teeth.off' - # left temporomandibular joint - 1 vertex (.ply/.off/coordinates) - p['axis_pt1'] = [0.1458893, -73.13895, 227.3909] - # right temporomandibular joint - 1 vertex (.ply/.off/coordinates) - p['axis_pt2'] = f'{path}/RTMJ.off' p['muscles'] = [ { @@ -43,12 +38,30 @@ def parms(d={}): 'method': 'T' } ] - p['food'] = 'fixteeth' # type of food: 'fixteeth' or 'cylinder' - p['fixations'] = { - 'contact_pts': ['x', 'y', 'z'], # constrained degrees of freedom - 'axis_pt1': ['x', 'y', 'z'], - 'axis_pt2': ['x', 'y', 'z'] - } + + # list of fixed nodes + # with "name" a description string used for output purposes, + # "nodes" a list of nodes (coordinates) or a mesh file, + # "direction" a list of constrained degrees of freedom + + p['fixations'] = [ + { + 'name': 'contact_pts', + 'nodes': [[-10.20362, -17.46838, -229.9061], + [-11.92466, 26.3042, -229.5354]], + 'direction': 'x' + }, + { + 'name': 'axis_pt1', + 'nodes': f'{path}/LTMJ.off', + 'direction': ['x', 'y', 'z'] + }, + { + 'name': 'axis_pt2', + 'nodes': [-8.716309, 79.13171, 233.8385], + 'direction': ['x', 'z'] + } + ] # material properties p['Young'] = 17000. # [MPa] elastic modulus - bone: 17-20 GPa @@ -91,10 +104,6 @@ def solve(p={}): nods_pos = np.reshape(nods_pos, (len(nods_no), 3)) # print(f'nods_pos={nods_pos}') - # load other surface meshes => groups['axis_pt1'], ... - for meshname in ['axis_pt1', 'axis_pt2']: - create_group(p[meshname], nods_no, nods_pos, groups, meshname) - # load muscle groups mgroups = {} # stores muscle group data and loads for muscle in p['muscles']: @@ -121,10 +130,9 @@ def solve(p={}): 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, - nods_pos, groups, 'contact_pts') + # create groups of fixed nodes + for fix in p['fixations']: + create_group(fix['nodes'], nods_no, nods_pos, groups, fix['name']) # -------------------------------------------------------------------------- @@ -136,14 +144,14 @@ def solve(p={}): # material do_not_delete.append( fem.Medium(pbl, "all", E=p['Young'], nu=p['Poisson']) ) + # add Dirichlet conditions (fixations) # print('== ADDING DIRICHLETs...') - # axis of rotation - for gname in ['axis_pt1', 'axis_pt2']: - for d in p['fixations'][gname]: - do_not_delete.append(fem.Dirichlet(pbl, gname, dof=d.upper(), val=0.0) ) + for fix in p['fixations']: + for d in fix['direction']: + do_not_delete.append(fem.Dirichlet(pbl, fix['name'], dof=d.upper(), val=0.0) ) + # add Neumann conditions (loads) from muscle loads # print('== ADDING NEUMANNs...') - # mshpts = geometry.getMesh().getPointSet() for name, gr in mgroups.items(): # print(f'len(gr.ntags)={len(gr.ntags)}') for i,load in enumerate(gr.loads): @@ -151,20 +159,9 @@ def solve(p={}): do_not_delete.append(fem.Neumann(pbl, int(gr.ntags[i]), 'Y', load.x[1]) ) # TODO: regler "dof=" qui marche pas do_not_delete.append(fem.Neumann(pbl, int(gr.ntags[i]), 'Z', load.x[2]) ) - # type of food - - if p['food'] == 'fixteeth': - # teeth are fixed along chosen directions - for d in p['fixations']['contact_pts']: - do_not_delete.append(fem.Dirichlet(pbl, 'contact_pts', dof=d.upper(), val=0.0) ) - else: - raise Exception('food = "fixteeth" only supported') - - # python only! # do_not_delete.append(fem.Dirichlet(pbl, "all", dof='T', val=0.0) ) - # Time integration scheme # print (pbl) # print('SOLVE.......') diff --git a/models/dolicorhynchops/dolicorhynchops_10k.py b/models/dolicorhynchops/dolicorhynchops_10k.py index 5a3ef24..99e195c 100644 --- a/models/dolicorhynchops/dolicorhynchops_10k.py +++ b/models/dolicorhynchops/dolicorhynchops_10k.py @@ -9,10 +9,6 @@ def parms(d={}): import os path = os.path.join(os.path.dirname(__file__),'10k') p['bone'] = f'{path}/mandible.stl' - p['contact_pts'] = [[-10.20362, -17.46838, -229.9061], - [-11.92466, 26.3042, -229.5354]] - p['axis_pt1'] = f'{path}/LTMJ.off' - p['axis_pt2'] = [-8.716309, 79.13171, 233.8385] p['muscles'] = [ { 'file': f'{path}/Lmuscle.stl', @@ -28,11 +24,24 @@ def parms(d={}): 'method': 'T' } ] - p['fixations'] = { - 'contact_pts': ['x'], - 'axis_pt1': ['x', 'y', 'z'], - 'axis_pt2': ['x', 'z'] - } + p['fixations'] = [ + { + 'name': 'contact_pts', + 'nodes': [[-10.20362, -17.46838, -229.9061], + [-11.92466, 26.3042, -229.5354]], + 'direction': 'x' + }, + { + 'name': 'axis_pt1', + 'nodes': f'{path}/LTMJ.off', + 'direction': ['x', 'y', 'z'] + }, + { + 'name': 'axis_pt2', + 'nodes': [-8.716309, 79.13171, 233.8385], + 'direction': ['x', 'z'] + } + ] # material properties p['density'] = 1.850e-9 # [T/mm³] -- GitLab