From 9ddd2d8d87f658e251c9470252da3b588aabc0d9 Mon Sep 17 00:00:00 2001
From: Sam Calisch <s.calisch@gmail.com>
Date: Sun, 3 Dec 2017 18:04:01 -0500
Subject: [PATCH] finish adding torques
---
flexure/CMD | 2 +
flexure/flexure_stiffness.py | 85 ++++++++++++++++++++++++++++++------
2 files changed, 73 insertions(+), 14 deletions(-)
diff --git a/flexure/CMD b/flexure/CMD
index 25049fe..9d58a8b 100644
--- a/flexure/CMD
+++ b/flexure/CMD
@@ -1 +1,3 @@
python flexure_stiffness.py -Q -M simulate -f 15 -flexure_type cyclic -w .0007 -t .0015 -l .010
+python flexure_stiffness.py -Q -M simulate -f 15 -flexure_type cyclic -w .0007 -t .0015 -l .010
+python flexure_stiffness.py -Q -M simulate -f 15 -flexure_type cyclic -w .0006 -t .0025 -l .009 -bd 10
diff --git a/flexure/flexure_stiffness.py b/flexure/flexure_stiffness.py
index ba1910f..90870a6 100644
--- a/flexure/flexure_stiffness.py
+++ b/flexure/flexure_stiffness.py
@@ -20,6 +20,16 @@ def plot_connections(nodes,beamsets):
ax.plot(nodes[seg,0], nodes[seg,1], nodes[seg,2], c=cmap(i))
plt.show()
+
+def get_rotation_from_nodes(nodes,axis,disps,targets):
+ base = asarray(axis[0]);
+ a = axis[1]-base;
+ a = a/magnitudes(a) #set up axis coordinates
+ x = nodes[targets]-base
+ d = x - dot(x, a)[...,None]*a
+ y = [dot(a,b) for a,b in zip(disps[targets,:3], cross(d,a))]
+ return arctan2(y,magnitudes(d)),magnitudes(d)
+
def run_frame3dd(args,nodes,global_args,beam_sets,constraints,loads):
write_frame3dd_file(nodes,global_args,beam_sets,constraints,loads)
cmd = ["frame3dd", "-i",global_args['frame3dd_filename']+'.csv']
@@ -64,6 +74,15 @@ def build(args):
])
if args.flexure_type == 'cyclic':
beams = vstack((beams, array([[4,20],[9,21],[14,22],[19,23]])))
+
+ nodes = vstack((nodes, array([ [args.sensor_radius,0,0],[0,args.sensor_radius,0],[-args.sensor_radius,0,0],[0,-args.sensor_radius,0] ])))
+ solid_beams = vstack(( solid_beams, array([
+ [24,0],[24,10],[24,20],[24,22],
+ [25,0],[25,10],[25,21],[25,23],
+ [26,5],[26,15],[26,21],[26,23],
+ [27,5],[27,15],[27,20],[27,22]
+ ])))
+
return nodes, beams, solid_beams
def run_simulation(args):
@@ -71,12 +90,13 @@ def run_simulation(args):
nodes,beams,solid_beams = build(args)
global_args = {
'n_modes':args.n_modes,'length_scaling':args.length_scaling,'exagerration':10,
- 'node_radius':zeros(shape(nodes)[0]),'frame3dd_filename':args.base_filename+"_frame3dd"
+ 'zoom_scale':2.,'node_radius':zeros(shape(nodes)[0]),
+ 'frame3dd_filename':args.base_filename+"_frame3dd"
}
clean_up_frame3dd(global_args['frame3dd_filename'])
beam_sets = [
(beams,{'E':args.E,'nu':args.nu,'rho':args.rho,'cross_section':'rectangular','d2':args.w,'d1':args.t,'roll':0.,'loads':[],'beam_divisions':args.bd,'prestresses':[]}),
- (solid_beams,{'E':10*args.E,'nu':args.nu,'rho':args.rho,'cross_section':'rectangular','d1':args.w,'d2':args.t,'roll':0.,'loads':[],'beam_divisions':args.bd,'prestresses':[]})
+ (solid_beams,{'E':10*args.E,'nu':args.nu,'rho':args.rho,'cross_section':'rectangular','d1':.003,'d2':.003,'roll':0.,'loads':[],'beam_divisions':1,'prestresses':[]})
]
if args.flexure_type == 'mirrored':
fixed_nodes = [2,4,7,9,12,14,17,19]
@@ -85,24 +105,58 @@ def run_simulation(args):
constraints = [{'node':node,'DOF':dof,'value':0} for dof in [0,1,2,3,4,5] for node in fixed_nodes]
+ loaded_nodes = [0,5,10,15,20,21,22,23]
+ sensor_nodes = [24,25,26,27]
+ results = []
+
for force_dof in [0,1,2]:
- loaded_nodes = [0,5,10,15]
loads = [{'node':n,'DOF':force_dof,'value':args.force/len(loaded_nodes)} for n in loaded_nodes]
-
run_frame3dd(args,nodes,global_args,beam_sets,constraints,loads)
- results = {}
- results['beam_mass'] = compute_mass(nodes,beam_sets)
disps = read_frame3dd_displacements(global_args['frame3dd_filename'])
force_disp = average(disps[loaded_nodes,force_dof])
- print "Degree of freedom: %d"%force_dof
- print "Force applied: %.1f N"%(args.force)
- print "Calculated displacement: %.1f microns"%(force_disp*1e6)
-
- #todo: read average displacement of loaded nodes
- #todo: plot displacements vs. design parameters
+ #print "Degree of freedom: %d"%force_dof
+ #print "Force applied: %.1f N"%(args.force)
+ #print "Displacement at sensor: %.1f microns"%(force_disp*1e6)
+ results.append( {'dof':force_dof, 'force/torque':args.force, 'displacement':force_disp} )
+ #sys.exit(0)
+ #Tx, Ty
+ torque_force = args.torque/(.5*args.sep)/len(loaded_nodes) # F = torque / dist
+ for torque_dof in [0,1]:
+ loads = [{'node':n,'DOF':torque_dof,'value':torque_force if nodes[n][2]>0 else -torque_force} for n in loaded_nodes]
+ run_frame3dd(args,nodes,global_args,beam_sets,constraints,loads)
+ disps = read_frame3dd_displacements(global_args['frame3dd_filename'])
+ moving_sensor_nodes = [24,26] if torque_dof==0 else [25,27]
+
+ #print disps[sensor_nodes]
+ #axis = (array([0,0,0]), array([0,-1,0]) if torque_dof==0 else array([1,0,0]) )
+ #rots,ds = get_rotation_from_nodes(nodes,axis,disps/global_args['length_scaling'],loaded_nodes)
+ #print average(rots)
+ #print "Degree of freedom: %d"%torque_dof
+ #print "Torque applied: %.1f Nm"%(args.torque)
+ #print "Torque force: %.2f N"%torque_force
+ #print "Displacement at sensor: %.1f microns"%(sin(average(rots))*args.sensor_radius*1e6)
+ #print "Displacement at sensor: %.1f microns"%( average(abs(disps[moving_sensor_nodes,2]))*1e6 )
+ results.append( {'dof':torque_dof+3, 'force/torque':args.torque, 'displacement':average(abs(disps[moving_sensor_nodes,2]))} )
+
+ #Tz
+ torque_force = args.torque/args.attach_radius/len(loaded_nodes) # F = torque / dist
+ def node_to_force(n):
+ return array([-nodes[n,1],nodes[n,0],0]) / sqrt(nodes[n,0]**2 + nodes[n,1]**2)
+ loads = [{'node':n,'DOF':i,'value':torque_force*node_to_force(n)[i]} for n in loaded_nodes for i in [0,1]]
+ run_frame3dd(args,nodes,global_args,beam_sets,constraints,loads)
+ disps = read_frame3dd_displacements(global_args['frame3dd_filename'])
+ #axis = (array([0,0,0]), array([0,0,1]))
+ #rots,ds = get_rotation_from_nodes(nodes,axis,disps/global_args['length_scaling'],loaded_nodes)
+ #print average(rots)
+ #print "Degree of freedom: %d"%2
+ #print "Torque applied: %.1f Nm"%(args.torque)
+ #print "Torque force: %.2f N"%torque_force
+ #print "Displacement at sensor: %.1f microns"%(sin(average(rots))*args.sensor_radius*1e6)
+ #print "Displacement at sensor: %.1f microns"%( average(magnitudes(disps[sensor_nodes,:3]))*1e6 )
+ results.append( {'dof':5, 'force/torque':args.torque, 'displacement':average(magnitudes(disps[sensor_nodes,:3]))} )
- #results['fundamental_frequency'] = read_lowest_mode(global_args['frame3dd_filename']+'.csv')
+ #todo: plot displacements vs. design parameters
return results
def find_stability_threshold(args):
@@ -146,7 +200,8 @@ if __name__ == '__main__':
parser.add_argument('-M','--mode',choices=('simulate','search', 'visualize'), required=True)
parser.add_argument('-flexure_type','--flexure_type',choices=('cyclic','mirrored'), required=True)
parser.add_argument('-Q','--quiet',action='store_true',help='Whether to suppress frame3dd output')
- parser.add_argument("-f","--force", type=double, default=.1, help="force to apply")
+ parser.add_argument("-f","--force", type=double, default=.1, help="force to apply (N)")
+ parser.add_argument("-torque","--torque", type=double, default=1., help="torque to apply (Nm)")
#parser.add_argument("-fr","--force_res", type=double, default=.01, help="Final resolution of force for search mode")
parser.add_argument("-w","--w", type=double, default=.0005, help="width of flexure (m)")
@@ -154,6 +209,7 @@ if __name__ == '__main__':
parser.add_argument("-l","--l", type=double, default=.0068, help="length of flexure segment (m)")
parser.add_argument("-attach_radius","--attach_radius", type=double, default=.0043, help="distance from z axis to flexure attachment (m)")
parser.add_argument("-sep","--sep", type=double, default=.025, help="flexure plate z separation (m)")
+ parser.add_argument("-sensor_radius","--sensor_radius", type=double, default=.012, help="distance from rotation axis to sensor (m)")
parser.add_argument("-bd","--bd", type=int, default=1, help='how many divisions for each rod, useful in buckling analysis')
parser.add_argument("-E","--E", type=double, default=70e9, help="Young's Modulus of laminate")
@@ -171,6 +227,7 @@ if __name__ == '__main__':
print "Critical stress: %.3f MPa"%(last_res['stress']/1e6)
elif args.mode=='simulate':
res = run_simulation(args)
+ print res
#print "Fundamental frequency: %.3f Hz"%res['fundamental_frequency']
#print "Stress: %.3f MPa"%(res['stress']/1e6)
elif args.mode=='visualize':
--
GitLab