Robust (and probabilistic) design (MF2024)
Lecture 2
MF2024@md.kth.se http://www.kth.se/itm/inst/mmk/edu/inst_kurser/md/MF2024/
1
Today
Monte Carlo Simulation Introduction to Individual Assignment 1 Introduction to L3a and L3b Screw Tightening Exercise
Reliability
- Machine reliability Costs of unexpected stand still in production lines are extremely high (steel mill line 10.000/hour, car paint shop plant 100.000/hour). - Product reliability Data storage devices
Reliability engineering
Reliability engineering
Component reliability System reliability
Procedures to establish component reliability require probability analysis.
Procedures to establish / optimize system reliability are
Fault tree analysis (FTA) Failure Mode and Effect Analysis (FMEA) Root Cause Failure analysis
Component Reliability
ISO: If a joint is fixed using a nut and bolt from the same grade there is no need to size the nut while the standard length of the nut is such sized that the screw will fail before the thread is stripped.
Case: Thread stripping strength
Component Reliability
Case: Thread stripping strength
Rule of thumb: Ath= 0.5d0L0 , d0 =(d2+d3)/2 Standard nut height 0.8d: L0= 0.75d M8: d2=7.188 mm, d3=6.466 mm, L0=6 mm Ath= 0.5d0L0 = 64.3 mm
2
Component Reliability
Case: Thread stripping strength
ISO: If a joint is fixed using a nut and bolt from the same grade there is no need to size the nut while the standard length of the nut is such sized that the screw will fail before the thread is stripped. The shear strength of a material is often expressed in the ratio of the shear to tensile strength, for ductile materials like steel /=0.580. M8: At = (/4)d0 , Ath = (1/0.58)At = 63.1 mm Rule of thumb: Ath= 0.5d0L0 = 64.3 mm
Rule of thumb corresponds with ISO
2 2 2
Component Reliability
Case: Thread stripping strength
Calculated shear strength: Rm = 360 MPa = 0.58Rm = 208.8 MPa L0 = 4 mm, M8 F = (As/0.58)(4/0.75d)=8.8kN
ISO
Measured shear strength: F = 18 kN
If the shear strength is much higher than needed it might be possible to select a thinner sheet resulting in lower cost and weight.
Component Reliability
Classification of failures Class 1 2 3 4 5 6 7 8 Range [kN] 11.50 - 13.25 13.25 - 15.00 15.00 - 16.75 16.75 18.50 18.50 20.25 20.25 22.00 22.00 23.75 23.75 25.50 ni 1 3 6 13 18 13 6 1 61
Case: Thread stripping strength
10
Component Reliability
Classification of failures Class Range [kN] 1 2 3 4 5 6 7 8 11.50 - 13.25 13.25 - 15.00 15.00 - 16.75 16.75 18.50 18.50 20.25 20.25 22.00 22.00 23.75 23.75 25.50 ni 1 3 6 13 18 13 6 1 f(t)=ni/n 0.016 0.049 0.098 0.213 0.295 0.213 0.098 0.016
Case: Thread stripping strength
61 1 The probability of failure F(t<=15kN)=1/61+3/61=0.016+0.049=0.0656 (6.6%)
11
Component Reliability
Summary
Case: Thread stripping strength Always check the validity of Rules of thumb Material properties may vary The deterministic approach might be useful with a large factor of safety. High tech reliability designs require a probabilistic approach
12
Normal distribution tolerance interval.
13
Component Reliability
Case: Hydraulic cylinder
Containerized Ring Cranes
Lifting capacity 1600 MT Own weight 600 MT Length main mast 97 m Height approximate 165 m Maximum ballast 1500 MT
14
Component Reliability
Case: Hydraulic cylinder
15
Component Reliability
Case: Hydraulic cylinder
16
Component Reliability
Case: Hydraulic cylinder dpiston = 299.5 300.0, p = 299.75, p = 0.5/6 dgroove = 294.5 295.0, g = 294.75, g = 0.5/6 hring = 2.9 3.1, r = 3.0, r = 0.2/6 Radial clearance Cr:
Worst case Cr = dgroove,min/2 + hring,min - dpiston,max/2 = 0.15 mm Cr = g/2 + r - p/2 = 0.5 mm Assumed ring deformation = 0.1 mm Measured ring deformation = 0.3 mm
17
Evolution of Design Technology
Trial & Error Probabilistic Empirical Mathematical
Deterministic (Factors of Safety)
Stochastic (Risk Quantified)
Random Experimentation Experience-based
Graphical Approaches Systematic Experimentation
Models based on system physics Point estimates
18
Simulations based on system physics Robust Solutions
18
19
Typical objective
20
Main methods for predicting response variation
21
What is Monte-Carlo Simulation?
An experiment performed on a computer rather than performed in an engineering laboratory If a system parameter (e.g. a design parameter) is known to follow a certain probability distribution, the performance of the system is studied by considering several possible values of the parameter Do not require normal distributed parameters
22
Work flow
23
Matlab
% % % % % % % % % % % % % % % % % % % % n % % % % Example Monte Carlo Simulation in Matlab Hydraulic cylinder case 1.1 in van Beek page 14-15 radial clearence rc is determined by a determistic and a probabalistic calculation rc=dgroove/2+hring-dpiston/2 rcmax=maximum allowable radial clearance uniform distribution case reference van Beek Advanced Engineering Design rc=radial clearence (mm) dgroove=diameter groove (mm) hring=thickness ring(mm) dpiston=diameter piston(mm) Generate n samples from a normal distribution r = ( rand(n,1) * sd ) + mu mu : mean sd : standard deviation Generate n samples from a uniform distribution r = a + rand(n,1) * (b-a) a : minimum b : maximum = 100000; % The number of function evaluations --- Generate vectors of random inputs dgroove ~ Uniform distribution U(a=294.5 b=295.0) hring ~ Uniform distribution U(a=2.9 b=3.1) dpiston ~ Uniform distribution U(a=299.5 b=300.0)
24
dgroove = 294.5 + rand(n,1) * ( 295.0 - 294.5 ); hring = 2.9 + rand(n,1) * ( 3.1 - 2.9 ); dpiston = 299.5 + rand(n,1) * ( 300.0 - 299.5 ); dgroovem=294.75; hringm=3.0; dpistonm=299.75; % --- Run the simulation crm=dgroovem/2+hringm-dpistonm/2 crmin=min(dgroove)/2+min(hring)-max(dpiston)/2 cr=dgroove/2+hring-dpiston/2; % --- Create a histogram of the results (50 bins) figure(1) hist(dgroove,50); title('diameter groove (mm)') figure(2) hist(hring,50); title('thickness ring (mm)') figure(3) hist(dpiston,50); title('diameter piston (mm)') figure(4) hist(cr,50); title('radial clearence (mm)') % --- Calculate summary statistics cr_mean = mean(cr) cr_std = std(cr) cr_median = median(cr) %Count fraction of points below cr_median fracmedian = sum(cr <= cr_median)/n %Count fraction of points below lower acceptable radial clearance 0.3 mm fraclow = sum(cr <= 0.3)/n
25
Simulation results
diameter piston (mm) 2500 2000 1500
1000
500
0 299.5 299.55 299.6 299.65 299.7 299.75 299.8 299.85 299.9 299.95
300
diameter groove (mm) 2500 2500
thickness ring (mm)
2000
2000
1500
1500
1000
1000
500
500
0 294.5 294.55 294.6 294.65 294.7 294.75 294.8 294.85 294.9 294.95
295
0 2.85
2.9
2.95
3.05
3.1
26
Simulation results
radial clearence (mm) 4500 4000 3500 3000 2500 2000 1500 1000 500 0 0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
fraction of points below lower acceptable radial clearance 0.3 mm = 4.5 %
27
Mating force in snap fastener
Common examples of integrated fasteners. (a) Module with four cantilever lugs; (b) cover with two cantilever and two rigid lugs; (c) separable snap joint for chassis cover.
28
Mating force in snap fastener
bh 2 E max P= 6L
tan + W =P 1 tan
3EI P= 3 L
29
Stochastic input variables
Normal distributed secant modulus of elasticity
m = 2.2 GPa = 0.1 GPa
Uniform distributed coefficient of friction
min = 0.54 max = 0.66
How is the mating force affected?
30
Tasks
Use Monte Carlo simulation to calculate the distribution of the mating force. Calculate the normal distribution parameters for the mating force (mean and standard deviation). Calculate the rejection rate given the maximum allowable mating force 60 N.
31
Useful matlab commands: rand
n-by-m matrix containing pseudorandom values drawn from the standard uniform distribution on the open interval (0,1): r = rand(n,m); n-by-m matrix containing pseudorandom values drawn from the standard uniform distribution on the open interval (min,max): r = min + (max-min)*rand(n,m);
32
Useful matlab commands: randn
n-by-m matrix containing pseudorandom values drawn from the standard normal distribution: r = randn(n,m); n-by-m matrix containing pseudorandom values drawn from the normal distribution with mean value mu and standard deviation sigma: r = mu + sigma*randn(n,m);
33
Monte Carlo simulations using Matlab
% Normal distributed secant modulus of elasticity [MPa] Emean = 2.2e3; Estd = 0.1e3; % Parameters % Uniform distributed coefficient of friction [-] mumin = 0.54; mumax = 0.66; % Parameters N = 100; % Number of trials for i = 1:N % Generate pseudorandom instances for input variables E(i) = Emean + Estd*randn; mu(i) = mumin + (mumax-mumin)*rand; % Calculate instances for output variables P(i) = b*h^2*E(i)*epsilon/(6*L); W(i) = P(i)*(tan(alpha)+mu(i))/(1-mu(i)*tan(alpha)); end
34
Alternative code
Use element-by-element operations (.* and ./) instead of for-loop
% Normal distributed secant modulus of elasticity [MPa] Emean = 2.2e3; Estd = 0.1e3; % Parameters % Uniform distributed coefficient of friction [-] mumin = 0.54; mumax = 0.66; % Parameters % Generate array with random instances for input variables N = 100; % Number of trials mu = mumin + (mumax-mumin)*rand(N,1); E = Emean + Estd*randn(N,1); % Calculate arrays with instances for output variables P = b*h^2*E*epsilon/(6*L); W = P.*(tan(alpha)+mu)./(1-mu*tan(alpha));
35
Trials
3000 0.8 0.6 2000
E [GPa]
[-]
1000 0
0.4 0.2 0
50 Trial
100
150
50 Trial
100
150
80 60
W [N]
40 20 0
50 Trial
100
150
36
Histogram
20 15 10 5 0 1.8 20 15 10 5 0 0.5
2.2 E [GPa]
2.4
2.6
0.55
0.6 [-]
0.65
0.7
20 15 10 5 0 45
50
55 W [N]
60
65
37
Normal distribution
Probability Distribution Function 0.2 0.15
pdf [-]
0.1 0.05 0 30
35
40
55 60 65 W [N] Cumulative Distribution Function
45
50
70
75
80
cdf [-]
0.5
0 30
35
40
45
50
55 W [N]
60
65
70
75
80
38
Rejection rate
Criteria: Reject if W>60 N
Probability Distribution Function 0.2 0.15
Probability of acceptance Probability of rejection
pdf [-]
0.1 0.05 0 30
35
40
45
55 60 65 W [N] Cumulative Distribution Function
50
70
75
80
0,90
cdf [-]
Rejection rate: (1-0.90) = 0.10
0.5
0 30
35
40
45
50
55 W [N]
60
65
70
75
80
39
Coefficients of friction for common snap fastener polymers, reference Bayer snap fit joints for plastics, A design guide, Bayer Polymers Division, 1966.
40
Use Monte Carlo simulation to calculate the distribution of the mating force and the rejection rate given the maximum allowable mating force 60 N. Also calculate the normal distribution parameters for the mating force (mean and standard deviation). Furthermore you should suggest alternative materials that can meet a demand of rejection rate below 1%.
Figure. Cantilever snap joint, figure Anders Sderberg KTH Machine Design.
41
Matlab
% Example Monte Carlo Simulation in Matlab % Snap fastener mating force for a cantilever snap joint with rectangular % cross section and constant thickness made of polycarbonate %W=P*(my+tan(alfa))/(1-my*tan(alfa)) % P=b*h^2*E*epsi/(6*L) % reference Hamrock, Jacobson and Schmid Fundamentals of machine elements % %W=mating force (N) %P=deformation force (N) %my=coefficient of friction (-) %alfa=wedge angle of cantilever tip (rad) %b=cantilever length (m) %h=cantilever thickness (%) %E=secant modulus of elasticity %epsi=allowable maximum strain for polymer material %L=cantilever length (m) % % % % % % % % % % Generate n samples from a normal distribution r = ( randn(n,1) * sd ) + mu mu : mean sd : standard deviation Generate n samples from a uniform distribution r = a + rand(n,1) * (b-a) a : minimum b : maximum
42
Matlab
n = 100000; % The number of function evaluations % --- Generate vectors of random inputs % E ~ Normal distribution N(mean=2.2e9,sd=0.1e9) % CES edupack E = 2 - 2.4 GPa % my ~ Uniform distribution U(a=0.54,b=0.66) %Self mated polymer Byaer Snap-fits for plastics smaller friction coefficent E = ( randn(n,1) * 0.1e9 ) + 2.2e9; my = 0.54 + rand(n,1) * ( 0.66 - 0.54 ); b=0.0095; L=0.019; h=0.0035; alfa=20*pi/180; epsi=0.02; % --- Run the simulation % Note the use of element-wise multiplication P = b*h^2.*E*epsi/6/L; W = P.*(my+tan(alfa))./(1-my*tan(alfa));
43
Matlab
% --- Create a histogram of the results (50 bins) figure(1) hist(E,50); title('secant modulus of elasticity (Pa)') figure(2) hist(my,50); title('coefficient of friction') figure(3) hist(P,50); title('deformation force (N)') figure(4) hist(W,50); title('mating force (N)') % --- Calculate summary statistics W_mean = mean(W) W_std = std(W) %Count fraction of points above upper acceptable mating force 60 N frachigh = sum(W >= 60)/n
44
Results fraction of points above upper acceptable mating force 60 N: 13 %, mean=55.4 N and std=3.8 N
mating force (N) 6000 5000
4000
3000
2000
1000
0 40
45
50
55
60
65
70
75
45
Ways to improve!
Choose material with lower friction coefficient Choose a material with lower Young's modulus
46
Coefficients of friction for common snap fastener polymers, reference Bayer snap fit joints for plastics, A design guide, Bayer Polymers Division, 1966.
47
Result polystrene
% --- Generate vectors of random inputs %polystyrene % E ~ Normal distribution N(mean=1.68e9,sd=0.1e9) % CES edupack E = 1.2 - 2.6 GPa % my ~ Uniform distribution U(a=0.48,b=0.60) %Self mated polymer Byaer Snap-fits for plastics smaller friction coefficent E = ( randn(n,1) * 0.1e9 ) + 1.68e9; my = 0.48 + rand(n,1) * ( 0.60 - 0.48 ); b=0.0095; L=0.019; h=0.0035; alfa=20*pi/180; epsi=0.02;
48
Result polystyrene
mating force (N) 6000 5000
4000
3000
2000
1000
0 25
30
35
40
45
50
55
49
Maximum allowable strain 1.6 % (polystyrene = 1.6 %) CES Edupack
% --- Generate vectors of random inputs %polystyrene % E ~ Normal distribution N(mean=1.68e9,sd=0.1e9) % CES edupack E = 1.2 - 2.6 GPa % my ~ Uniform distribution U(a=0.48,b=0.60) % maximum allovable strain 1.6 % CES EDUPACK %Self mated polymer Byaer Snap-fits for plastics smaller friction coefficent E = ( randn(n,1) * 0.1e9 ) + 1.68e9; my = 0.48 + rand(n,1) * ( 0.60 - 0.48 ); b=0.0095; L=0.019; h=0.0035; alfa=20*pi/180; epsi=0.016;
50
Change angle to 15 degree
% --- Generate vectors of random inputs % E ~ Normal distribution N(mean=2.2e9,sd=0.1e9) % CES edupack E = 2 - 2.4 GPa % my ~ Uniform distribution U(a=0.54,b=0.66) %Self mated polymer Byaer Snap-fits for plastics smaller friction coefficent E = ( randn(n,1) * 0.1e9 ) + 2.2e9; my = 0.54 + rand(n,1) * ( 0.66 - 0.54 ); b=0.0095; L=0.019; h=0.0035; alfa=15*pi/180; epsi=0.02;
51
Angle 15 degree zero rejection rate
mating force (N) 6000 5000
4000
3000
2000
1000
0 30
35
40
45
50
55
60
52
Assignment 1: Threaded fastener
-a probabilistic study in MF2024,
2011
53
Threaded fastener
The tightening torque required to obtained an specified axial force depends on: Friction in thread Friction in head contact
Fmin Fi Fmax
54
Tasks
Monte Carlo simulation of the distribution of the axial force and an analysis of the rejection rate, due to variations in tightening torque and coefficients of friction Correlation of model with experimental measurements to find the variations in the coefficients of friction (i.e. inverse modeling) Modify the tightening torque in the simulation model to reduce the rejection rate.
55
About the assignment
Detailed instruction on course homepage Requires experimental results (L3 Screw tightening exercise) This assignment is graded and should be solved individually Present the result in a written technical report Send the report to mf2024@md.kth.se no later than Friday November 18, 2011
56
Grading of individual assignments
Each individual assignment gives 0-20 points, with a minimum of 10 points required to pass The written exam gives 0-40 points, with 15 points required to pass Final grading is based on the sum of the points for the individual assignments and the written examination
57
Laboratory experiment (L3)
Measurement of axial force in threaded fastener Thursday 3/11 & Monday7/11 4 groups of 10 students 3/11 (8.15-9, 9.15-10), 7/11 (10.15-11, 11.15-12) Sign up on the list that is passed around Machine Design, Lab Hall, Brinellvgen 83, ground floor
58
What is next?
L3 Screw tightening exercise - Group 1, 2 Thursday 3/11 - Group 3, 4 Monday 7/11 L4 Supervision Assignment 1 - Thursday 10/11
59
