1% MBDyn (C) is a multibody analysis code.
2% http://www.mbdyn.org
3%
4% Copyright (C) 1996-2013
5%
6% Pierangelo Masarati	<masarati@aero.polimi.it>
7% Paolo Mantegazza	<mantegazza@aero.polimi.it>
8%
9% Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
10% via La Masa, 34 - 20156 Milano, Italy
11% http://www.aero.polimi.it
12%
13% Changing this copyright notice is forbidden.
14%
15% This program is free software; you can redistribute it and/or modify
16% it under the terms of the GNU General Public License as published by
17% the Free Software Foundation (version 2 of the License).
18%
19%
20% This program is distributed in the hope that it will be useful,
21% but WITHOUT ANY WARRANTY; without even the implied warranty of
22% MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
23% GNU General Public License for more details.
24%
25% You should have received a copy of the GNU General Public License
26% along with this program; if not, write to the Free Software
27% Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
28
29% AUTHOR: Reinhard Resch <r.resch@secop.com>
30%        Copyright (C) 2011(-2013) all rights reserved.
31%
32%        The copyright of this code is transferred
33%        to Pierangelo Masarati and Paolo Mantegazza
34%        for use in the software MBDyn as described
35%        in the GNU Public License version 2.1
36
37function elem = BallBearingContact(pMbElem, pDM, HP)
38    elem.pMbElem = pMbElem;
39    elem.pDM = pDM;
40
41	if ( ~HP.IsKeyWord("ballradius") )
42		error("ball bearing contact(%d): keyword \"ball radius\" expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
43	endif
44
45	elem.R = HP.GetReal();
46
47	if ( HP.IsKeyWord("elasticmodulus") )
48		elem.E = HP.GetReal();
49	else
50		if ( ~HP.IsKeyWord("elasticmodulusball") )
51			error("ball bearing contact(%d): keyword \"elastic modulus ball\" expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
52		endif
53
54		E1 = HP.GetReal();
55
56		if ( ~HP.IsKeyWord("poissonratioball") )
57			error("ball bearing contact(%d): keyword \"poisson ratio ball\" expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
58		endif
59
60		nu1 = HP.GetReal();
61
62		if ( ~HP.IsKeyWord("elasticmoduluswasherdisk") )
63			error("ball bearing contact(%d): keyword \"elastic modulus washer disk\" expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
64		endif
65
66		E2 = HP.GetReal();
67
68		if ( ~HP.IsKeyWord("poissonratiowasherdisk") )
69			error("ball bearing contact(%d): keyword \"poisson ratio washer disk\" expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
70		endif
71
72		nu2 = HP.GetReal();
73
74		elem.E = 1 / ((1 - nu1^2) / E1 + (1 - nu2^2) / E2);
75	endif
76
77    if ( ~HP.IsKeyWord("node1") )
78        error("ball bearing contact(%d): keyword node1 expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
79    endif
80
81    elem.pNode1 = pDM.ReadNode(HP, "STRUCTURAL");
82
83    if ( ~HP.IsKeyWord("node2") )
84        error("ball bearing contact(%d): keyword node2 expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
85    endif
86
87    elem.pNode2 = pDM.ReadNode(HP, "STRUCTURAL");
88
89    if ( HP.IsKeyWord("offset") )
90        elem.o2 = HP.GetPosRel(elem.pNode2);
91    else
92        elem.o2 = zeros(3,1);
93    endif
94
95    if ( HP.IsKeyWord("orientation") )
96        elem.Rt2 = HP.GetRotRel(elem.pNode2);
97    else
98        elem.Rt2 = eye(3);
99    endif
100
101    if ( ~(HP.IsKeyWord("coulombfrictioncoefficient") || HP.IsKeyWord("coulombfrictioncoefficientx")) )
102        error("ball bearing contact(%d): keyword coulomb friction coefficient x expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
103    endif
104
105    mukx = HP.GetReal();
106
107    if ( mukx <= 0 )
108        error("ball bearing contact(%d): coulomb friction coefficient x must be greater than zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
109    endif
110
111    if ( HP.IsKeyWord("coulombfrictioncoefficienty") )
112        muky = HP.GetReal();
113    else
114        muky = mukx;
115    endif
116
117    if ( muky <= 0 )
118        error("ball bearing contact(%d): coulomb friction coefficient y must be greater than zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
119    endif
120
121    elem.Mk = diag([mukx, muky]);
122
123    if ( ~(HP.IsKeyWord("staticfrictioncoefficient") || HP.IsKeyWord("staticfrictioncoefficientx")) )
124        error("ball bearing contact(%d): keyword static friction coefficient x expected at line %d", pMbElem.GetLabel(), HP.GetLineData());
125    endif
126
127    musx = HP.GetReal();
128
129    if ( musx <= 0 )
130        error("ball bearing contact(%d): static friction coefficient x must be greater than zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
131    endif
132
133    if ( HP.IsKeyWord("staticfrictioncoefficienty") )
134        musy = HP.GetReal();
135    else
136        musy = musx;
137    endif
138
139    if ( musy <= 0 )
140        error("ball bearing contact(%d): static friction coefficient y must be greater than zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
141    endif
142
143    elem.Ms = diag([musx, musy]);
144
145    if ( ~HP.IsKeyWord("slidingvelocitycoefficient") )
146        error("ball bearing contact(%d): keyword sliding velocity coefficient expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
147    endif
148
149    elem.vs = HP.GetReal();
150
151    if ( elem.vs <= 0 )
152        error("ball bearing contact(%d): sliding velocity coefficient must be greater than zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
153    endif
154
155    if ( HP.IsKeyWord("slidingvelocityexponent") )
156        elem.gamma = HP.GetReal();
157    else
158        elem.gamma = 1;
159    endif
160
161    if ( ~(HP.IsKeyWord("microslipstiffness") || HP.IsKeyWord("microslipstiffnessx")) )
162        error("ball bearing contact(%d): keyword micro slip stiffness x expected at line %s", pMbElem.GetLabel(), HP.GetLineData());
163    endif
164
165    sigma0x = HP.GetReal();
166
167    if ( sigma0x <= 0 )
168        error("ball bearing contact(%d): micro slip stiffness x must be greater than zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
169    endif
170
171    if ( HP.IsKeyWord("microslipstiffnessy") )
172        sigma0y = HP.GetReal();
173    else
174        sigma0y = sigma0x;
175    endif
176
177    if ( sigma0y <= 0 )
178        error("ball bearing contact(%d): micro slip stiffness y must be greater than zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
179    endif
180
181    elem.sigma0 = diag([sigma0x, sigma0y]);
182
183    if ( HP.IsKeyWord("microslipdamping") || HP.IsKeyWord("microslipdampingx") )
184        sigma1x = HP.GetReal();
185
186        if ( sigma1x < 0 )
187            error("ball bearing contact(%d): micro slip damping x must be greater than or equal to zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
188        endif
189
190        if ( HP.IsKeyWord("microslipdampingy") )
191            sigma1y = HP.GetReal();
192
193            if ( sigma1y < 0 )
194                error("ball bearing contact(%d): micro slip damping y must be greater than or equal to zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
195            endif
196        else
197            sigma1y = sigma1x;
198        endif
199    else
200        sigma1x = sigma1y = 0;
201    endif
202
203    elem.sigma1 = diag([sigma1x, sigma1y]);
204
205    if ( HP.IsKeyWord("viscousfrictioncoefficient") || HP.IsKeyWord("viscousfrictioncoefficientx") )
206        sigma2x = HP.GetReal();
207
208        if ( sigma2x < 0 )
209            error("ball bearing contact(%d): viscous friction coefficient x must be greater than or equal to zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
210        endif
211
212        if ( HP.IsKeyWord("viscousfrictioncoefficienty") )
213            sigma2y = HP.GetReal();
214
215            if ( sigma2y < 0 )
216                error("ball bearing contact(%d): viscous friction coefficient y must be greater than or equal to zero at line %s", pMbElem.GetLabel(), HP.GetLineData());
217            endif
218        else
219            sigma2y = sigma2x;
220        endif
221    else
222        sigma2x = sigma2y = 0;
223    endif
224
225    elem.sigma2 = diag([sigma2x, sigma2y]);
226
227    elem.z = zeros(2,1);
228    elem.zP = zeros(2,1);
229
230	if (HP.IsKeyWord("initialstictionstate"))
231		for i=1:2
232			elem.z(i) = HP.GetReal();
233		endfor
234	endif
235
236	if (HP.IsKeyWord("initialstictionstatederivative"))
237		for i=1:2
238			elem.zP(i) = HP.GetReal();
239		endfor
240	endif
241
242    elem.rolling_friction = "no";
243
244    if (HP.IsKeyWord("rollingfrictionmodel"))
245        if (HP.IsKeyWord("SKF"))
246            elem.rolling_friction = "SKF";
247            if (~HP.IsKeyWord("staticloadcapacity"))
248                error("ball bearing contact(%d): keyword \"staticloadcapacity\" expected at line %s\n", pMbElem.GetLabel(), HP.GetLineData());
249            endif
250            elem.C0 = HP.GetReal();
251        elseif (HP.IsKeyWord("rigidball"))
252            elem.rolling_friction = "rigid ball";
253        elseif (HP.IsKeyWord("no"))
254            elem.rolling_friction = "no";
255        else
256            error("ball bearing contact(%d): friction model not supported at line %s\n", pMbElem.GetLabel(), HP.GetLineData());
257        endif
258    endif
259
260    elem = class(elem, "BallBearingContact");
261endfunction
262