beef.general

General purpose functions.

  1'''
  2General purpose functions.
  3'''
  4
  5import numpy as np
  6from scipy.linalg import null_space as null
  7
  8
  9def n2d_ix(node_ixs, n_dofs=6):
 10    r'''
 11    Convert node indices to dof indices.
 12
 13    Arguments
 14    ----------
 15    node_ixs : int
 16        list/array of indices of nodes
 17    n_dofs : 6, optional    
 18        number of dofs per node
 19
 20    Returns
 21    ----------
 22    dof_ixs : int
 23        array of DOF indices corresponding to input node indices
 24    '''
 25    return np.hstack([node_ix*n_dofs+np.arange(n_dofs) for node_ix in node_ixs])
 26
 27
 28def extract_dofs(mat, dof_ix=[0,1,2], n_dofs=6, axis=0):
 29    r'''
 30    Extract selected dofs from matrix.
 31
 32    Arguments
 33    ------------
 34    mat : float
 35        numpy array (matrix) to extract DOFs from
 36    dof_ix : int, optional
 37        list of DOF indices (local numbering of DOFs in node) to extract
 38    n_dofs : 6, optional
 39        number of DOFs per node
 40    axis : 0, optional
 41        axis to pick along
 42
 43    Returns
 44    ----------
 45    mat_sub : float
 46        subselection of matrix
 47    '''
 48    get_dofs = np.vstack([np.arange(dof, mat.shape[axis],n_dofs) for dof in dof_ix]).T.flatten()
 49    return mat.take(get_dofs, axis=axis)
 50
 51def convert_dofs(dofs_in, node_type='beam3d', sort_output=True):
 52    r'''
 53    Convert string DOFs to indices.
 54
 55    Arguments
 56    ----------
 57    dofs_in : {'all', 'trans', 'rot'} or int
 58        string or list of indices
 59    node_type : {'beam2d', 'beam3d'}, optional
 60        type of nodes (defining number of DOFs per node)
 61    sort_output : True, optional
 62        whether or not to sort the DOF indices output
 63    
 64    Returns
 65    --------
 66    dofs_out : int
 67        array of DOF indices 
 68    '''
 69    dof_dict = dict(beam2d=dict(all=np.arange(0,3), trans=np.arange(0,2), rot=np.arange(2,3)), 
 70                        beam3d=dict(all=np.arange(0,6), trans=np.arange(0,3), rot=np.arange(3,6)))
 71
 72    if dofs_in in ['all', 'trans', 'rot']:
 73        dofs_out = dof_dict[node_type][dofs_in]
 74    else:
 75        dofs_out = dofs_in
 76        
 77    if sort_output:
 78        dofs_out = np.sort(dofs_out)
 79        
 80    return dofs_out
 81        
 82
 83def convert_dofs_list(dofs_list_in, n_nodes, node_type='beam3d', sort_output=True):
 84    r'''
 85    Convert DOFs list to correct format.
 86
 87    Arguments
 88    -----------
 89    dofs_list_in : str,int
 90        list of strings and indices corresponding to DOFs to
 91    n_nodes : int
 92        number of nodes the DOFs are related to
 93    node_type : {'beam3d', 'beam2d'}
 94        type of nodes (defining number of DOFs per node)
 95    sort_output : True
 96        whether or not to sort output
 97
 98    Returns
 99    -----------
100    dofs_list_out : int
101        array of DOF indices
102    '''
103    
104    contains_strings = np.any([type(d) is str for d in dofs_list_in])
105    
106    if type(dofs_list_in) is not list:  # single input (all, rot or trans)
107        dofs_list_out = [dofs_list_in]*n_nodes      
108    elif ~contains_strings and (len(dofs_list_in)<=6) and (np.max(dofs_list_in)<6):   #
109        dofs_list_out = [dofs_list_in]*n_nodes
110    elif len(dofs_list_in)!=n_nodes:
111        raise TypeError('Wrong input format of "dofs"')
112    else:
113        dofs_list_out = dofs_list_in
114        
115    for ix, dofs_list_out_i in enumerate(dofs_list_out):
116        dofs_list_out[ix] = convert_dofs(dofs_list_out_i, node_type=node_type, sort_output=sort_output)
117    
118    return dofs_list_out
119
120
121def transform_unit(e1, e2=None, e3=None, warnings=False):
122    r'''
123    Establish transformation matrix from e1 and temporary e2 or e3 vectors.
124
125    Arguments
126    -----------
127    e1 : float
128        unit vector describing element longitudinal direction
129    e2 : float, optional
130        temporary unit vector describing a chosen vector that's perpendicular to the longitudinal direction (approximate y-direction)
131    e3 : float, optional
132        temporary unit vector describing a chosen vector that's perpendicular to the longitudinal direction (approximate z-direction)
133        if both e2 and e3 are different from None, e2 is used (e3 disregarded)
134
135    Returns
136    -----------
137    T : float
138        transformation matrix
139
140
141    Notes
142    -----------
143    Input vectors \(\{e_1\}\) and \(\{e_{2p}\}\) are used to establish a third unit vector, as the cross product of the two, as follows:
144
145    $$
146    \{e_3\} = \{e_1\} \\times \{e_{2p}\}
147    $$
148
149    Then, a final vector is created based on a second cross product:
150    $$
151    \{e_2\} = \{e_3\} \\times \{e_1\}
152    $$
153
154    Finally, the transformation matrix is established as follows:
155
156    $$
157    [T] = \\begin{bmatrix}
158        \{e_1\}^T \\\\
159        \{e_2\}^T \\\\
160        \{e_3\}^T
161    \\end{bmatrix}
162    $$
163
164    where the unit vectors established are the rows of the transformation matrix.
165
166    '''
167
168    e1 = np.array(e1).flatten()
169
170    if (e2 is not None) and (e3 is not None):
171        e3 = None
172    
173    if e2 is not None:
174        e2 = np.array(e2).flatten()
175        e3_tmp = np.cross(e1, e2)         # Direction of the third unit vector
176        if np.all(e3==0):
177            if warnings:
178                print('Warning: e1 and e2 identical. Check orientations carefully!')
179
180            e2 = e2 + 0.1
181            e3 = np.cross(e1, e2) 
182        
183        e2 = np.cross(e3_tmp, e1)          # Direction of the second unit vector
184
185        e1 = e1/np.linalg.norm(e1)     # Normalize the direction vectors to become unit vectors
186        e2 = e2/np.linalg.norm(e2)
187        e3 = np.cross(e1, e2)
188
189    elif e3 is not None:
190        e3 = np.array(e3).flatten()
191        e2_tmp = np.cross(e3, e1)         # Direction of the third unit vector
192        e3 = np.cross(e1, e2_tmp)
193        e1 = e1/np.linalg.norm(e1)     # Normalize the direction vectors to become unit vectors
194        e3 = e3/np.linalg.norm(e3)
195        e2 = np.cross(e3, e1) 
196    
197    if e2 is None and e3 is None:
198        raise ValueError('Specify either e2 or e3')
199
200    T = np.vstack([e1,e2,e3])
201    
202    return T
203
204
205def gdof_ix_from_nodelabels(all_node_labels, node_labels, dof_ix=[0,1,2]):     # from nlfe - not debugged
206    '''
207    Get global DOF indices from node labels.
208
209    Arguments
210    ------------
211    all_node_labels : int
212        list or array of all node labels
213    node_labels : int
214        list or array of the node labels to get the indices of
215    dof_ix : int
216        local DOF indices to use
217
218    Returns
219    ---------
220    gdof_ix : int
221        global indices of DOFs of requested nodes
222    '''
223    if type(node_labels) is not list:
224        node_labels = [node_labels]
225    
226    node_ix = [np.where(nl==all_node_labels)[0] for nl in node_labels]
227    gdof_ix = gdof_from_nodedof(node_ix, dof_ix)
228    
229    return gdof_ix
230
231
232def gdof_from_nodedof(node_ix, dof_ixs, n_dofs=3, merge=True):
233    r'''
234    Get global DOF from node DOF.
235
236    Arguments
237    ----------
238    node_ix : int
239        node indices to establish global DOF indices of
240    dof_ixs : int
241        local DOF indices to include
242    n_dofs : 3, optional
243        number of DOFs (all)
244    merge : True, optional
245        whether or not to merge the global DOF indices from all nodes, to a single list
246
247    Returns
248    ---------
249    gdof_ix : int
250        global indices of DOFs of requested nodes
251
252    '''
253    gdof_ix = []
254    
255    if type(node_ix) is not list:
256        node_ix = [node_ix]
257    
258    if type(dof_ixs) is not list:
259        dof_ixs = [dof_ixs]*len(node_ix)
260        
261    elif len(dof_ixs) != len(node_ix):
262        dof_ixs = [dof_ixs]*len(node_ix)
263
264    for ix, n in enumerate(node_ix):
265        gdof_ix.append(n*n_dofs + np.array(dof_ixs[ix]))
266    
267    if merge:
268        gdof_ix = np.array(gdof_ix).flatten()
269        
270    return gdof_ix
271
272def B_to_dofpairs(B, master_val=1):
273    r'''
274    Establish pairs of indices of DOFs to couple from compatibility matrix B.
275
276    Arguments
277    ----------
278    B : float   
279        Lagrange compatiblity \([B]\) matrix
280    master_val : float
281        value used to identify master DOF
282
283    Returns
284    ----------
285    dof_pairs : int
286        indices of DOFs paired through input \([B]\) matrix
287
288    '''
289    n_constr, n_dofs = B.shape
290    dof_pairs = [None]*n_constr
291    
292    for icon in range(0, n_constr):
293        master = np.where(B[icon,:] == master_val)[0][0]
294        slave = np.where(B[icon,:] == -master_val)[0]
295        if len(slave) != 0:
296            slave = slave[0]
297        else:
298            slave = None
299        
300        dof_pairs[icon] = [master, slave]
301        dof_pairs = np.array(dof_pairs).T
302
303    return dof_pairs
304
305
306def dof_pairs_to_Linv(dof_pairs, n_dofs):    
307    r'''
308    Establish quasi-inverse of L from dof pairs describing constraints.
309
310    Arguments
311    ----------
312    dof_pairs : int   
313        list of lists of indices of DOFs paired
314    n_dofs : int
315        number of DOFs
316
317    Returns
318    ----------
319    Linv : int
320
321    Notes
322    ----------
323    \([L_{inv}]\) is a constructed quantity, conventient as it can construct the constrainted u (unique, free DOFs) from the full as follows:
324
325    $$
326    \{u_c\} = [L_{inv}] \{u\}
327    $$
328
329    It is worth noting that it only works for fully hard-constrained systems (no weak connections, or similar).
330
331    '''
332
333    slave_nodes = dof_pairs[dof_pairs[:,1]!=None, 1]
334    grounded_nodes = dof_pairs[dof_pairs[:,1]==None, 0]
335    
336    remove = np.unique(np.hstack([grounded_nodes, slave_nodes]))   
337    all_dofs = np.arange(0,n_dofs)
338    primal_dofs = np.setdiff1d(all_dofs, remove)
339    
340    Linv = np.zeros([len(primal_dofs), n_dofs])
341
342    for primal_dof_ix, primal_dof in enumerate(primal_dofs):
343        Linv[primal_dof_ix, primal_dof] = 1
344    
345    return Linv
346
347
348def compatibility_matrix(dof_pairs, n_dofs):
349    r'''
350    Establish compatibility matrix from specified pairs of DOFs.
351
352    Arguments
353    ----------
354    dof_pairs : int   
355        list of lists of indices of DOFs paired
356    n_dofs : int
357        number of DOFs in full system; defines the number of rows of \([B]\)
358
359    Returns
360    ----------
361    B : int
362        numpy array describing compatibility matrix \([B]\)
363    '''
364    n_constraints = dof_pairs.shape[0]
365    B = np.zeros([n_constraints, n_dofs])
366
367    for constraint_ix, dof_pair in enumerate(dof_pairs):
368        if dof_pair[1] is None:    # connected dof is None --> fixed to ground
369            B[constraint_ix, dof_pair[0]] = 1
370        else:
371            B[constraint_ix, dof_pair[0]] = 1
372            B[constraint_ix, dof_pair[1]] = -1
373    return B
374
375
376def lagrange_constrain(mat, dof_pairs, null=False):
377    r'''
378    Lagrange constrain matrix from specified DOF pairs.
379
380    Arguments
381    -----------------
382    mat : float
383        matrix to constrain
384    dof_pairs : int
385        list of lists of indices of DOFs paired
386    null : False, optional
387        to create a 0 matrix of the correct size, this can be set to True
388
389    Returns
390    ----------
391    mat_fixed : float
392        constrained matrix
393    '''
394    # Lagrange multiplier method - expand matrix
395    n_constraints = len(dof_pairs)
396    n_dofs = mat.shape[0]
397    
398    if not null:
399        B = compatibility_matrix(dof_pairs, n_dofs)
400    else:
401        B = np.zeros([n_constraints, n_dofs])
402        
403    O = np.zeros([n_constraints, n_constraints])
404    mat_fixed = np.vstack([np.hstack([mat, B.T]), np.hstack([B, O])])
405   
406    return mat_fixed
407
408def lagrange_constrain_from_B(mat, B, null=False):
409    r'''
410    Lagrange constrain matrix from specified compatibility matrix.
411
412    Arguments
413    -----------------
414    mat : float
415        matrix to constrain
416    B : int
417        numpy array describing compatibility matrix \([B]\)
418    null : False, optional
419        to create a 0 matrix of the correct size, this can be set to True
420
421    Returns
422    ----------
423    mat_fixed : float
424        constrained matrix
425    '''
426    # Lagrange multiplier method - expand matrix
427    n_constraints = B.shape[0]        
428    O = np.zeros([n_constraints, n_constraints])
429    mat_fixed = np.vstack([np.hstack([mat, B.T]), np.hstack([B, O])])
430   
431    return mat_fixed
432
433
434def expand_vector_with_zeros(vec, n_constraints):
435    r'''
436    Append vector with zeros based on number of constraints (n_constraints).
437    '''
438    vec_expanded = np.hstack([vec.flatten(), np.zeros(n_constraints)])[np.newaxis, :].T
439    return vec_expanded
440
441
442def blkdiag(mat, n):
443    return np.kron(np.eye(n), mat)
444
445
446def nodes_to_beam_element_matrix(node_labels, first_element_label=1):
447    r'''
448    Establish element matrix definition assuming ordered node labels.
449
450    Arguments
451    -----------
452    node_labels : int
453        list of integer node labels
454    first_element_label : int
455        first integer used in element matrix
456
457    Returns
458    -----------
459    element_matrix : int
460        element definition where each row represents an element as [element_label_i, node_label_1, node_label_2]
461
462    '''
463    n_nodes = len(node_labels)
464    n_elements = n_nodes-1
465    
466    element_matrix = np.empty([n_elements, 3])
467    element_matrix[:, 0] = np.arange(first_element_label,first_element_label+n_elements)
468    element_matrix[:, 1] = node_labels[0:-1]
469    element_matrix[:, 2] = node_labels[1:]
470
471    return element_matrix
472
473
474def get_phys_modes(phi, B, lambd=None, return_as_ix=False):
475    r'''
476    Get physical modes.
477
478    Arguments
479    -----------
480    phi : float
481    B : int
482        numpy array describing compatibility matrix \([B]\)
483    lambd : float, optional
484        standard value None does not 
485    return_as_ix : False, optional
486        whether or not to output as indices instead of eigenvectors and eigenvalues
487
488    Returns
489    ----------
490    lambd_phys : float
491        physical eigenvalue array
492    phi_phys : float
493        physical eigenvector array (modal transformation matrix)
494    phys_ix : int
495        indices of selected DOFs (if `return_as_ix=True` is passed)
496    '''
497    n_lagr = B.shape[0]
498    u_ = phi[0:-n_lagr, :]     #physical dofs
499    l_ = phi[-n_lagr:,:]       #lagr dofs
500    
501    # Compatibility at interface
502    mask1 = np.all((B @ u_)==0, axis=0)
503    
504    # Equilibrium at interface
505    L = null(B)
506    g = -B.T @ l_
507    mask2 = np.all((L.T @ g)==0, axis=0)
508    
509    phys_ix = np.logical_and(mask1, mask1)
510    
511    if return_as_ix:
512        return phys_ix
513    else:  
514        if lambd is None:
515            lambd_phys = None
516        else:
517            lambd_phys = lambd[phys_ix]
518            
519        phi_phys = u_[:, phys_ix]
520        
521        return lambd_phys, phi_phys    
522    
523    
524def ensure_list(list_in, levels=1, increase_only=True):
525    r'''
526    Ensure input variable is list.
527
528    Arguments
529    ----------
530    list_in : float, int
531        list or scalar
532    levels : int
533        number of levels of list (nest-level)
534    increase_only : True, optional
535        if True, only increase amount of list-wrapping
536
537    Returns
538    ----------
539    list_out : float, int
540        list
541    '''
542
543    dimensions = np.array(list_in).ndim 
544
545    if type(list_in) is not list:
546        list_out = [list_in]
547        dimensions = dimensions+1
548    else:
549        list_out = list_in * 1
550    
551    
552    if not increase_only:
553        while dimensions>levels:
554            list_out = list_out[0]
555            dimensions = dimensions - 1
556
557    while dimensions<levels:            
558        list_out = [list_out]
559        dimensions = dimensions + 1
560
561    return list_out
562
563
564def feature_matrix(master_dofs, values, slave_dofs=None, ndofs=None, return_small=False):
565    r'''
566    Arguments
567    ---------
568    master_dofs : int
569        list of indices of master DOFs
570    values : float
571        list of amplitudes/values for features
572    slave_dofs : int, optional
573        list of indices of slave DOFs (standard value None, makes grounded connections from master_dofs)
574    ndofs : int, optional
575        number of modes used to construct full matrix (not relevant if `return_small=True`)
576    return_small : False, optional
577        whether or not to return small or full matrix - if `return_small=True`, returns small ndofs-by-ndofs matrix and corresponding indices in global matrix; if `return_small=False` returns big ndofs_global-by-ndofs_global matrix
578    
579    Returns
580    ---------
581    mat : float
582        feature matrix
583    dof_ixs : int
584        list of indices of DOFs for relevant feature connectivity (only returned if `return_small=True`)
585
586    '''
587    
588    ndofs_small = len(master_dofs)
589    if type(values) is float or type(values) is int:
590        values = [float(values)]*ndofs_small
591    elif len(values) != len(master_dofs):
592        raise ValueError('Length of master_dofs and values must match')
593
594    if slave_dofs is None:
595        small = np.diag(values)
596        slave_dofs = []
597    else:
598        if len(slave_dofs) != len(master_dofs):
599            raise ValueError('Length of master_dofs and slave_dofs must match')
600        small = np.diag(np.hstack([values,values]))
601        for ix in range(ndofs_small):
602            small[ix, ix+ndofs_small] = small[ix+ndofs_small, ix] = -values[ix]
603
604    dof_ixs = np.hstack([master_dofs, slave_dofs]).astype(int)
605
606    if return_small:
607        return small, dof_ixs
608    else:
609        if len(set(dof_ixs)) != len(dof_ixs):
610            raise ValueError('Non-unique dof indices are provided')
611        if ndofs is None:
612            ndofs = np.max(dof_ixs)+1
613
614        large = np.zeros([ndofs, ndofs])
615        large[np.ix_(dof_ixs, dof_ixs)] = small
616        
617        return large
618
619def basic_coupled():
620    return np.array([[1, -1], [-1, 1]])
621
622def generic_beam_mat(L, yy=0, yz=0, yt=0, zy=0, zz=0, zt=0, ty=0, tz=0, tt=0):
623    mat = np.zeros([12,12])
624
625    mat[0:6, 0:6] = np.array([
626        [0,         0,          0,          0,          0,              0           ],
627        [0,         156*yy,    156*yz,    147*yt,    -22*L*yz,      22*L*yy    ],
628        [0,         156*zy,    156*zz,    147*zt,    -22*L*zz,      22*L*zy    ],
629        [0,         147*ty,    147*tz,    140*tt,     -21*L*tz,      21*L*ty   ],
630        [0,         -22*L*zy,  -22*L*zz,  -21*L*zt,  4*L**2*zz,     -4*L**2*zy ],
631        [0,         22*L*yy,   22*L*yz,   21*L*yt,   -4*L**2*yz,    4*L**2*yy  ],
632    ])
633
634    mat[0:6, 6:12] = np.array([
635        [0,         0,          0,          0,          0,              0            ],
636        [0,         54*yy,    54*yz,      63*yt,     13*L*yz,       -13*L*yy    ],
637        [0,         54*zy,    54*zz,      63*zt,     13*L*zz,       -13*L*zy    ],
638        [0,         63*ty,    63*tz,      70*tt,     14*L*tz,       -14*L*ty    ],
639        [0,         -13*L*zy,  -13*L*zz,  -14*L*zt,  -3*L**2*zz,     3*L**2*zy  ],
640        [0,         13*L*yy,   13*L*yz,   14*L*yt,   3*L**2*yz,     -3*L**2*yy  ],
641    ])
642
643    mat[6:12, 0:6] = np.array([
644        [0,         0,          0,          0,          0,              0            ],
645        [0,         54*yy,    54*yz,      63*yt,     -13*L*yz,       13*L*yy    ],
646        [0,         54*zy,    54*zz,      63*zt,     -13*L*zz,       13*L*zy    ],
647        [0,         63*ty,    63*tz,      70*tt,     -14*L*tz,       14*L*ty    ],
648        [0,         13*L*zy,  13*L*zz,    14*L*zt,   -3*L**2*zz,     3*L**2*zy  ],
649        [0,         -13*L*yy, -13*L*yz,   -14*L*yt,   3*L**2*yz,     -3*L**2*yy ],
650    ])
651
652    mat[6:12,6:12] = np.array([
653        [0,         0,          0,          0,          0,              0               ],
654        [0,         156*yy,    156*yz,    147*yt,    22*L*yz,      -22*L*yy        ],
655        [0,         156*zy,    156*zz,    147*zt,    22*L*zz,      -22*L*zy        ],
656        [0,         147*ty,    147*tz,    140*tt,    21*L*tz,      -21*L*ty        ],
657        [0,         22*L*zy,   22*L*zz,   21*L*zt,    4*L**2*zz,   -4*L**2*zy      ],
658        [0,         -22*L*yy,   -22*L*yz,   -21*L*yt,   -4*L**2*yz,    4*L**2*yy   ],
659    ])
660
661    return mat * L / 420
662
663
664def bar_foundation_stiffness(L, kx, ky, kz):    #only z and y currently, will be extended!
665    mat = np.vstack([ 
666        [kx*1/4, 0, 0, 0, 0, 0,     kx*1/4, 0, 0, 0, 0, 0],     #ux1
667        [0, ky*1/3, 0, 0, 0, 0,     0, ky*1/6, 0, 0, 0, 0],     #uy1
668        [0, 0, kz*1/3, 0, 0, 0,     0, 0, kz*1/6, 0, 0, 0],     #uz1
669
670        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #rx1
671        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #ry1
672        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #rz1
673
674        [kx*1/4, 0, 0, 0, 0, 0,     kx*1/4, 0, 0, 0, 0, 0],     #ux2
675        [0, ky*1/6, 0, 0, 0, 0,     0, ky*1/3, 0, 0, 0, 0],     #uy2
676        [0, 0, kz*1/6, 0, 0, 0,     0, 0, kz*1/3, 0, 0, 0],     #uz2
677
678        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #rx2
679        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #ry2
680        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0]]) * L      #rz2
681    return mat
def n2d_ix(node_ixs, n_dofs=6):
10def n2d_ix(node_ixs, n_dofs=6):
11    r'''
12    Convert node indices to dof indices.
13
14    Arguments
15    ----------
16    node_ixs : int
17        list/array of indices of nodes
18    n_dofs : 6, optional    
19        number of dofs per node
20
21    Returns
22    ----------
23    dof_ixs : int
24        array of DOF indices corresponding to input node indices
25    '''
26    return np.hstack([node_ix*n_dofs+np.arange(n_dofs) for node_ix in node_ixs])

Convert node indices to dof indices.

Arguments

node_ixs : int list/array of indices of nodes n_dofs : 6, optional
number of dofs per node

Returns
  • dof_ixs (int): array of DOF indices corresponding to input node indices
def extract_dofs(mat, dof_ix=[0, 1, 2], n_dofs=6, axis=0):
29def extract_dofs(mat, dof_ix=[0,1,2], n_dofs=6, axis=0):
30    r'''
31    Extract selected dofs from matrix.
32
33    Arguments
34    ------------
35    mat : float
36        numpy array (matrix) to extract DOFs from
37    dof_ix : int, optional
38        list of DOF indices (local numbering of DOFs in node) to extract
39    n_dofs : 6, optional
40        number of DOFs per node
41    axis : 0, optional
42        axis to pick along
43
44    Returns
45    ----------
46    mat_sub : float
47        subselection of matrix
48    '''
49    get_dofs = np.vstack([np.arange(dof, mat.shape[axis],n_dofs) for dof in dof_ix]).T.flatten()
50    return mat.take(get_dofs, axis=axis)

Extract selected dofs from matrix.

Arguments

mat : float numpy array (matrix) to extract DOFs from dof_ix : int, optional list of DOF indices (local numbering of DOFs in node) to extract n_dofs : 6, optional number of DOFs per node axis : 0, optional axis to pick along

Returns
  • mat_sub (float): subselection of matrix
def convert_dofs(dofs_in, node_type='beam3d', sort_output=True):
52def convert_dofs(dofs_in, node_type='beam3d', sort_output=True):
53    r'''
54    Convert string DOFs to indices.
55
56    Arguments
57    ----------
58    dofs_in : {'all', 'trans', 'rot'} or int
59        string or list of indices
60    node_type : {'beam2d', 'beam3d'}, optional
61        type of nodes (defining number of DOFs per node)
62    sort_output : True, optional
63        whether or not to sort the DOF indices output
64    
65    Returns
66    --------
67    dofs_out : int
68        array of DOF indices 
69    '''
70    dof_dict = dict(beam2d=dict(all=np.arange(0,3), trans=np.arange(0,2), rot=np.arange(2,3)), 
71                        beam3d=dict(all=np.arange(0,6), trans=np.arange(0,3), rot=np.arange(3,6)))
72
73    if dofs_in in ['all', 'trans', 'rot']:
74        dofs_out = dof_dict[node_type][dofs_in]
75    else:
76        dofs_out = dofs_in
77        
78    if sort_output:
79        dofs_out = np.sort(dofs_out)
80        
81    return dofs_out

Convert string DOFs to indices.

Arguments

dofs_in : {'all', 'trans', 'rot'} or int string or list of indices node_type : {'beam2d', 'beam3d'}, optional type of nodes (defining number of DOFs per node) sort_output : True, optional whether or not to sort the DOF indices output

Returns
  • dofs_out (int): array of DOF indices
def convert_dofs_list(dofs_list_in, n_nodes, node_type='beam3d', sort_output=True):
 84def convert_dofs_list(dofs_list_in, n_nodes, node_type='beam3d', sort_output=True):
 85    r'''
 86    Convert DOFs list to correct format.
 87
 88    Arguments
 89    -----------
 90    dofs_list_in : str,int
 91        list of strings and indices corresponding to DOFs to
 92    n_nodes : int
 93        number of nodes the DOFs are related to
 94    node_type : {'beam3d', 'beam2d'}
 95        type of nodes (defining number of DOFs per node)
 96    sort_output : True
 97        whether or not to sort output
 98
 99    Returns
100    -----------
101    dofs_list_out : int
102        array of DOF indices
103    '''
104    
105    contains_strings = np.any([type(d) is str for d in dofs_list_in])
106    
107    if type(dofs_list_in) is not list:  # single input (all, rot or trans)
108        dofs_list_out = [dofs_list_in]*n_nodes      
109    elif ~contains_strings and (len(dofs_list_in)<=6) and (np.max(dofs_list_in)<6):   #
110        dofs_list_out = [dofs_list_in]*n_nodes
111    elif len(dofs_list_in)!=n_nodes:
112        raise TypeError('Wrong input format of "dofs"')
113    else:
114        dofs_list_out = dofs_list_in
115        
116    for ix, dofs_list_out_i in enumerate(dofs_list_out):
117        dofs_list_out[ix] = convert_dofs(dofs_list_out_i, node_type=node_type, sort_output=sort_output)
118    
119    return dofs_list_out

Convert DOFs list to correct format.

Arguments

dofs_list_in : str,int list of strings and indices corresponding to DOFs to n_nodes : int number of nodes the DOFs are related to node_type : {'beam3d', 'beam2d'} type of nodes (defining number of DOFs per node) sort_output : True whether or not to sort output

Returns
  • dofs_list_out (int): array of DOF indices
def transform_unit(e1, e2=None, e3=None, warnings=False):
122def transform_unit(e1, e2=None, e3=None, warnings=False):
123    r'''
124    Establish transformation matrix from e1 and temporary e2 or e3 vectors.
125
126    Arguments
127    -----------
128    e1 : float
129        unit vector describing element longitudinal direction
130    e2 : float, optional
131        temporary unit vector describing a chosen vector that's perpendicular to the longitudinal direction (approximate y-direction)
132    e3 : float, optional
133        temporary unit vector describing a chosen vector that's perpendicular to the longitudinal direction (approximate z-direction)
134        if both e2 and e3 are different from None, e2 is used (e3 disregarded)
135
136    Returns
137    -----------
138    T : float
139        transformation matrix
140
141
142    Notes
143    -----------
144    Input vectors \(\{e_1\}\) and \(\{e_{2p}\}\) are used to establish a third unit vector, as the cross product of the two, as follows:
145
146    $$
147    \{e_3\} = \{e_1\} \\times \{e_{2p}\}
148    $$
149
150    Then, a final vector is created based on a second cross product:
151    $$
152    \{e_2\} = \{e_3\} \\times \{e_1\}
153    $$
154
155    Finally, the transformation matrix is established as follows:
156
157    $$
158    [T] = \\begin{bmatrix}
159        \{e_1\}^T \\\\
160        \{e_2\}^T \\\\
161        \{e_3\}^T
162    \\end{bmatrix}
163    $$
164
165    where the unit vectors established are the rows of the transformation matrix.
166
167    '''
168
169    e1 = np.array(e1).flatten()
170
171    if (e2 is not None) and (e3 is not None):
172        e3 = None
173    
174    if e2 is not None:
175        e2 = np.array(e2).flatten()
176        e3_tmp = np.cross(e1, e2)         # Direction of the third unit vector
177        if np.all(e3==0):
178            if warnings:
179                print('Warning: e1 and e2 identical. Check orientations carefully!')
180
181            e2 = e2 + 0.1
182            e3 = np.cross(e1, e2) 
183        
184        e2 = np.cross(e3_tmp, e1)          # Direction of the second unit vector
185
186        e1 = e1/np.linalg.norm(e1)     # Normalize the direction vectors to become unit vectors
187        e2 = e2/np.linalg.norm(e2)
188        e3 = np.cross(e1, e2)
189
190    elif e3 is not None:
191        e3 = np.array(e3).flatten()
192        e2_tmp = np.cross(e3, e1)         # Direction of the third unit vector
193        e3 = np.cross(e1, e2_tmp)
194        e1 = e1/np.linalg.norm(e1)     # Normalize the direction vectors to become unit vectors
195        e3 = e3/np.linalg.norm(e3)
196        e2 = np.cross(e3, e1) 
197    
198    if e2 is None and e3 is None:
199        raise ValueError('Specify either e2 or e3')
200
201    T = np.vstack([e1,e2,e3])
202    
203    return T

Establish transformation matrix from e1 and temporary e2 or e3 vectors.

Arguments

e1 : float unit vector describing element longitudinal direction e2 : float, optional temporary unit vector describing a chosen vector that's perpendicular to the longitudinal direction (approximate y-direction) e3 : float, optional temporary unit vector describing a chosen vector that's perpendicular to the longitudinal direction (approximate z-direction) if both e2 and e3 are different from None, e2 is used (e3 disregarded)

Returns
  • T (float): transformation matrix
Notes

Input vectors ({e_1}) and ({e_{2p}}) are used to establish a third unit vector, as the cross product of the two, as follows:

$$ {e_3} = {e_1} \times {e_{2p}} $$

Then, a final vector is created based on a second cross product: $$ {e_2} = {e_3} \times {e_1} $$

Finally, the transformation matrix is established as follows:

$$ [T] = \begin{bmatrix} {e_1}^T \\ {e_2}^T \\ {e_3}^T \end{bmatrix} $$

where the unit vectors established are the rows of the transformation matrix.

def gdof_ix_from_nodelabels(all_node_labels, node_labels, dof_ix=[0, 1, 2]):
206def gdof_ix_from_nodelabels(all_node_labels, node_labels, dof_ix=[0,1,2]):     # from nlfe - not debugged
207    '''
208    Get global DOF indices from node labels.
209
210    Arguments
211    ------------
212    all_node_labels : int
213        list or array of all node labels
214    node_labels : int
215        list or array of the node labels to get the indices of
216    dof_ix : int
217        local DOF indices to use
218
219    Returns
220    ---------
221    gdof_ix : int
222        global indices of DOFs of requested nodes
223    '''
224    if type(node_labels) is not list:
225        node_labels = [node_labels]
226    
227    node_ix = [np.where(nl==all_node_labels)[0] for nl in node_labels]
228    gdof_ix = gdof_from_nodedof(node_ix, dof_ix)
229    
230    return gdof_ix

Get global DOF indices from node labels.

Arguments

all_node_labels : int list or array of all node labels node_labels : int list or array of the node labels to get the indices of dof_ix : int local DOF indices to use

Returns
  • gdof_ix (int): global indices of DOFs of requested nodes
def gdof_from_nodedof(node_ix, dof_ixs, n_dofs=3, merge=True):
233def gdof_from_nodedof(node_ix, dof_ixs, n_dofs=3, merge=True):
234    r'''
235    Get global DOF from node DOF.
236
237    Arguments
238    ----------
239    node_ix : int
240        node indices to establish global DOF indices of
241    dof_ixs : int
242        local DOF indices to include
243    n_dofs : 3, optional
244        number of DOFs (all)
245    merge : True, optional
246        whether or not to merge the global DOF indices from all nodes, to a single list
247
248    Returns
249    ---------
250    gdof_ix : int
251        global indices of DOFs of requested nodes
252
253    '''
254    gdof_ix = []
255    
256    if type(node_ix) is not list:
257        node_ix = [node_ix]
258    
259    if type(dof_ixs) is not list:
260        dof_ixs = [dof_ixs]*len(node_ix)
261        
262    elif len(dof_ixs) != len(node_ix):
263        dof_ixs = [dof_ixs]*len(node_ix)
264
265    for ix, n in enumerate(node_ix):
266        gdof_ix.append(n*n_dofs + np.array(dof_ixs[ix]))
267    
268    if merge:
269        gdof_ix = np.array(gdof_ix).flatten()
270        
271    return gdof_ix

Get global DOF from node DOF.

Arguments

node_ix : int node indices to establish global DOF indices of dof_ixs : int local DOF indices to include n_dofs : 3, optional number of DOFs (all) merge : True, optional whether or not to merge the global DOF indices from all nodes, to a single list

Returns
  • gdof_ix (int): global indices of DOFs of requested nodes
def B_to_dofpairs(B, master_val=1):
273def B_to_dofpairs(B, master_val=1):
274    r'''
275    Establish pairs of indices of DOFs to couple from compatibility matrix B.
276
277    Arguments
278    ----------
279    B : float   
280        Lagrange compatiblity \([B]\) matrix
281    master_val : float
282        value used to identify master DOF
283
284    Returns
285    ----------
286    dof_pairs : int
287        indices of DOFs paired through input \([B]\) matrix
288
289    '''
290    n_constr, n_dofs = B.shape
291    dof_pairs = [None]*n_constr
292    
293    for icon in range(0, n_constr):
294        master = np.where(B[icon,:] == master_val)[0][0]
295        slave = np.where(B[icon,:] == -master_val)[0]
296        if len(slave) != 0:
297            slave = slave[0]
298        else:
299            slave = None
300        
301        dof_pairs[icon] = [master, slave]
302        dof_pairs = np.array(dof_pairs).T
303
304    return dof_pairs

Establish pairs of indices of DOFs to couple from compatibility matrix B.

Arguments

B : float
Lagrange compatiblity ([B]) matrix master_val : float value used to identify master DOF

Returns
  • dof_pairs (int): indices of DOFs paired through input ([B]) matrix
def dof_pairs_to_Linv(dof_pairs, n_dofs):
307def dof_pairs_to_Linv(dof_pairs, n_dofs):    
308    r'''
309    Establish quasi-inverse of L from dof pairs describing constraints.
310
311    Arguments
312    ----------
313    dof_pairs : int   
314        list of lists of indices of DOFs paired
315    n_dofs : int
316        number of DOFs
317
318    Returns
319    ----------
320    Linv : int
321
322    Notes
323    ----------
324    \([L_{inv}]\) is a constructed quantity, conventient as it can construct the constrainted u (unique, free DOFs) from the full as follows:
325
326    $$
327    \{u_c\} = [L_{inv}] \{u\}
328    $$
329
330    It is worth noting that it only works for fully hard-constrained systems (no weak connections, or similar).
331
332    '''
333
334    slave_nodes = dof_pairs[dof_pairs[:,1]!=None, 1]
335    grounded_nodes = dof_pairs[dof_pairs[:,1]==None, 0]
336    
337    remove = np.unique(np.hstack([grounded_nodes, slave_nodes]))   
338    all_dofs = np.arange(0,n_dofs)
339    primal_dofs = np.setdiff1d(all_dofs, remove)
340    
341    Linv = np.zeros([len(primal_dofs), n_dofs])
342
343    for primal_dof_ix, primal_dof in enumerate(primal_dofs):
344        Linv[primal_dof_ix, primal_dof] = 1
345    
346    return Linv

Establish quasi-inverse of L from dof pairs describing constraints.

Arguments

dof_pairs : int
list of lists of indices of DOFs paired n_dofs : int number of DOFs

Returns
  • Linv (int):
Notes

([L_{inv}]) is a constructed quantity, conventient as it can construct the constrainted u (unique, free DOFs) from the full as follows:

$$ {u_c} = [L_{inv}] {u} $$

It is worth noting that it only works for fully hard-constrained systems (no weak connections, or similar).

def compatibility_matrix(dof_pairs, n_dofs):
349def compatibility_matrix(dof_pairs, n_dofs):
350    r'''
351    Establish compatibility matrix from specified pairs of DOFs.
352
353    Arguments
354    ----------
355    dof_pairs : int   
356        list of lists of indices of DOFs paired
357    n_dofs : int
358        number of DOFs in full system; defines the number of rows of \([B]\)
359
360    Returns
361    ----------
362    B : int
363        numpy array describing compatibility matrix \([B]\)
364    '''
365    n_constraints = dof_pairs.shape[0]
366    B = np.zeros([n_constraints, n_dofs])
367
368    for constraint_ix, dof_pair in enumerate(dof_pairs):
369        if dof_pair[1] is None:    # connected dof is None --> fixed to ground
370            B[constraint_ix, dof_pair[0]] = 1
371        else:
372            B[constraint_ix, dof_pair[0]] = 1
373            B[constraint_ix, dof_pair[1]] = -1
374    return B

Establish compatibility matrix from specified pairs of DOFs.

Arguments

dof_pairs : int
list of lists of indices of DOFs paired n_dofs : int number of DOFs in full system; defines the number of rows of ([B])

Returns
  • B (int): numpy array describing compatibility matrix ([B])
def lagrange_constrain(mat, dof_pairs, null=False):
377def lagrange_constrain(mat, dof_pairs, null=False):
378    r'''
379    Lagrange constrain matrix from specified DOF pairs.
380
381    Arguments
382    -----------------
383    mat : float
384        matrix to constrain
385    dof_pairs : int
386        list of lists of indices of DOFs paired
387    null : False, optional
388        to create a 0 matrix of the correct size, this can be set to True
389
390    Returns
391    ----------
392    mat_fixed : float
393        constrained matrix
394    '''
395    # Lagrange multiplier method - expand matrix
396    n_constraints = len(dof_pairs)
397    n_dofs = mat.shape[0]
398    
399    if not null:
400        B = compatibility_matrix(dof_pairs, n_dofs)
401    else:
402        B = np.zeros([n_constraints, n_dofs])
403        
404    O = np.zeros([n_constraints, n_constraints])
405    mat_fixed = np.vstack([np.hstack([mat, B.T]), np.hstack([B, O])])
406   
407    return mat_fixed

Lagrange constrain matrix from specified DOF pairs.

Arguments

mat : float matrix to constrain dof_pairs : int list of lists of indices of DOFs paired null : False, optional to create a 0 matrix of the correct size, this can be set to True

Returns
  • mat_fixed (float): constrained matrix
def lagrange_constrain_from_B(mat, B, null=False):
409def lagrange_constrain_from_B(mat, B, null=False):
410    r'''
411    Lagrange constrain matrix from specified compatibility matrix.
412
413    Arguments
414    -----------------
415    mat : float
416        matrix to constrain
417    B : int
418        numpy array describing compatibility matrix \([B]\)
419    null : False, optional
420        to create a 0 matrix of the correct size, this can be set to True
421
422    Returns
423    ----------
424    mat_fixed : float
425        constrained matrix
426    '''
427    # Lagrange multiplier method - expand matrix
428    n_constraints = B.shape[0]        
429    O = np.zeros([n_constraints, n_constraints])
430    mat_fixed = np.vstack([np.hstack([mat, B.T]), np.hstack([B, O])])
431   
432    return mat_fixed

Lagrange constrain matrix from specified compatibility matrix.

Arguments

mat : float matrix to constrain B : int numpy array describing compatibility matrix ([B]) null : False, optional to create a 0 matrix of the correct size, this can be set to True

Returns
  • mat_fixed (float): constrained matrix
def expand_vector_with_zeros(vec, n_constraints):
435def expand_vector_with_zeros(vec, n_constraints):
436    r'''
437    Append vector with zeros based on number of constraints (n_constraints).
438    '''
439    vec_expanded = np.hstack([vec.flatten(), np.zeros(n_constraints)])[np.newaxis, :].T
440    return vec_expanded

Append vector with zeros based on number of constraints (n_constraints).

def blkdiag(mat, n):
443def blkdiag(mat, n):
444    return np.kron(np.eye(n), mat)
def nodes_to_beam_element_matrix(node_labels, first_element_label=1):
447def nodes_to_beam_element_matrix(node_labels, first_element_label=1):
448    r'''
449    Establish element matrix definition assuming ordered node labels.
450
451    Arguments
452    -----------
453    node_labels : int
454        list of integer node labels
455    first_element_label : int
456        first integer used in element matrix
457
458    Returns
459    -----------
460    element_matrix : int
461        element definition where each row represents an element as [element_label_i, node_label_1, node_label_2]
462
463    '''
464    n_nodes = len(node_labels)
465    n_elements = n_nodes-1
466    
467    element_matrix = np.empty([n_elements, 3])
468    element_matrix[:, 0] = np.arange(first_element_label,first_element_label+n_elements)
469    element_matrix[:, 1] = node_labels[0:-1]
470    element_matrix[:, 2] = node_labels[1:]
471
472    return element_matrix

Establish element matrix definition assuming ordered node labels.

Arguments

node_labels : int list of integer node labels first_element_label : int first integer used in element matrix

Returns
  • element_matrix (int): element definition where each row represents an element as [element_label_i, node_label_1, node_label_2]
def get_phys_modes(phi, B, lambd=None, return_as_ix=False):
475def get_phys_modes(phi, B, lambd=None, return_as_ix=False):
476    r'''
477    Get physical modes.
478
479    Arguments
480    -----------
481    phi : float
482    B : int
483        numpy array describing compatibility matrix \([B]\)
484    lambd : float, optional
485        standard value None does not 
486    return_as_ix : False, optional
487        whether or not to output as indices instead of eigenvectors and eigenvalues
488
489    Returns
490    ----------
491    lambd_phys : float
492        physical eigenvalue array
493    phi_phys : float
494        physical eigenvector array (modal transformation matrix)
495    phys_ix : int
496        indices of selected DOFs (if `return_as_ix=True` is passed)
497    '''
498    n_lagr = B.shape[0]
499    u_ = phi[0:-n_lagr, :]     #physical dofs
500    l_ = phi[-n_lagr:,:]       #lagr dofs
501    
502    # Compatibility at interface
503    mask1 = np.all((B @ u_)==0, axis=0)
504    
505    # Equilibrium at interface
506    L = null(B)
507    g = -B.T @ l_
508    mask2 = np.all((L.T @ g)==0, axis=0)
509    
510    phys_ix = np.logical_and(mask1, mask1)
511    
512    if return_as_ix:
513        return phys_ix
514    else:  
515        if lambd is None:
516            lambd_phys = None
517        else:
518            lambd_phys = lambd[phys_ix]
519            
520        phi_phys = u_[:, phys_ix]
521        
522        return lambd_phys, phi_phys    

Get physical modes.

Arguments

phi : float B : int numpy array describing compatibility matrix ([B]) lambd : float, optional standard value None does not return_as_ix : False, optional whether or not to output as indices instead of eigenvectors and eigenvalues

Returns
  • lambd_phys (float): physical eigenvalue array
  • phi_phys (float): physical eigenvector array (modal transformation matrix)
  • phys_ix (int): indices of selected DOFs (if return_as_ix=True is passed)
def ensure_list(list_in, levels=1, increase_only=True):
525def ensure_list(list_in, levels=1, increase_only=True):
526    r'''
527    Ensure input variable is list.
528
529    Arguments
530    ----------
531    list_in : float, int
532        list or scalar
533    levels : int
534        number of levels of list (nest-level)
535    increase_only : True, optional
536        if True, only increase amount of list-wrapping
537
538    Returns
539    ----------
540    list_out : float, int
541        list
542    '''
543
544    dimensions = np.array(list_in).ndim 
545
546    if type(list_in) is not list:
547        list_out = [list_in]
548        dimensions = dimensions+1
549    else:
550        list_out = list_in * 1
551    
552    
553    if not increase_only:
554        while dimensions>levels:
555            list_out = list_out[0]
556            dimensions = dimensions - 1
557
558    while dimensions<levels:            
559        list_out = [list_out]
560        dimensions = dimensions + 1
561
562    return list_out

Ensure input variable is list.

Arguments

list_in : float, int list or scalar levels : int number of levels of list (nest-level) increase_only : True, optional if True, only increase amount of list-wrapping

Returns
  • list_out (float, int): list
def feature_matrix(master_dofs, values, slave_dofs=None, ndofs=None, return_small=False):
565def feature_matrix(master_dofs, values, slave_dofs=None, ndofs=None, return_small=False):
566    r'''
567    Arguments
568    ---------
569    master_dofs : int
570        list of indices of master DOFs
571    values : float
572        list of amplitudes/values for features
573    slave_dofs : int, optional
574        list of indices of slave DOFs (standard value None, makes grounded connections from master_dofs)
575    ndofs : int, optional
576        number of modes used to construct full matrix (not relevant if `return_small=True`)
577    return_small : False, optional
578        whether or not to return small or full matrix - if `return_small=True`, returns small ndofs-by-ndofs matrix and corresponding indices in global matrix; if `return_small=False` returns big ndofs_global-by-ndofs_global matrix
579    
580    Returns
581    ---------
582    mat : float
583        feature matrix
584    dof_ixs : int
585        list of indices of DOFs for relevant feature connectivity (only returned if `return_small=True`)
586
587    '''
588    
589    ndofs_small = len(master_dofs)
590    if type(values) is float or type(values) is int:
591        values = [float(values)]*ndofs_small
592    elif len(values) != len(master_dofs):
593        raise ValueError('Length of master_dofs and values must match')
594
595    if slave_dofs is None:
596        small = np.diag(values)
597        slave_dofs = []
598    else:
599        if len(slave_dofs) != len(master_dofs):
600            raise ValueError('Length of master_dofs and slave_dofs must match')
601        small = np.diag(np.hstack([values,values]))
602        for ix in range(ndofs_small):
603            small[ix, ix+ndofs_small] = small[ix+ndofs_small, ix] = -values[ix]
604
605    dof_ixs = np.hstack([master_dofs, slave_dofs]).astype(int)
606
607    if return_small:
608        return small, dof_ixs
609    else:
610        if len(set(dof_ixs)) != len(dof_ixs):
611            raise ValueError('Non-unique dof indices are provided')
612        if ndofs is None:
613            ndofs = np.max(dof_ixs)+1
614
615        large = np.zeros([ndofs, ndofs])
616        large[np.ix_(dof_ixs, dof_ixs)] = small
617        
618        return large
Arguments

master_dofs : int list of indices of master DOFs values : float list of amplitudes/values for features slave_dofs : int, optional list of indices of slave DOFs (standard value None, makes grounded connections from master_dofs) ndofs : int, optional number of modes used to construct full matrix (not relevant if return_small=True) return_small : False, optional whether or not to return small or full matrix - if return_small=True, returns small ndofs-by-ndofs matrix and corresponding indices in global matrix; if return_small=False returns big ndofs_global-by-ndofs_global matrix

Returns
  • mat (float): feature matrix
  • dof_ixs (int): list of indices of DOFs for relevant feature connectivity (only returned if return_small=True)
def basic_coupled():
620def basic_coupled():
621    return np.array([[1, -1], [-1, 1]])
def generic_beam_mat(L, yy=0, yz=0, yt=0, zy=0, zz=0, zt=0, ty=0, tz=0, tt=0):
623def generic_beam_mat(L, yy=0, yz=0, yt=0, zy=0, zz=0, zt=0, ty=0, tz=0, tt=0):
624    mat = np.zeros([12,12])
625
626    mat[0:6, 0:6] = np.array([
627        [0,         0,          0,          0,          0,              0           ],
628        [0,         156*yy,    156*yz,    147*yt,    -22*L*yz,      22*L*yy    ],
629        [0,         156*zy,    156*zz,    147*zt,    -22*L*zz,      22*L*zy    ],
630        [0,         147*ty,    147*tz,    140*tt,     -21*L*tz,      21*L*ty   ],
631        [0,         -22*L*zy,  -22*L*zz,  -21*L*zt,  4*L**2*zz,     -4*L**2*zy ],
632        [0,         22*L*yy,   22*L*yz,   21*L*yt,   -4*L**2*yz,    4*L**2*yy  ],
633    ])
634
635    mat[0:6, 6:12] = np.array([
636        [0,         0,          0,          0,          0,              0            ],
637        [0,         54*yy,    54*yz,      63*yt,     13*L*yz,       -13*L*yy    ],
638        [0,         54*zy,    54*zz,      63*zt,     13*L*zz,       -13*L*zy    ],
639        [0,         63*ty,    63*tz,      70*tt,     14*L*tz,       -14*L*ty    ],
640        [0,         -13*L*zy,  -13*L*zz,  -14*L*zt,  -3*L**2*zz,     3*L**2*zy  ],
641        [0,         13*L*yy,   13*L*yz,   14*L*yt,   3*L**2*yz,     -3*L**2*yy  ],
642    ])
643
644    mat[6:12, 0:6] = np.array([
645        [0,         0,          0,          0,          0,              0            ],
646        [0,         54*yy,    54*yz,      63*yt,     -13*L*yz,       13*L*yy    ],
647        [0,         54*zy,    54*zz,      63*zt,     -13*L*zz,       13*L*zy    ],
648        [0,         63*ty,    63*tz,      70*tt,     -14*L*tz,       14*L*ty    ],
649        [0,         13*L*zy,  13*L*zz,    14*L*zt,   -3*L**2*zz,     3*L**2*zy  ],
650        [0,         -13*L*yy, -13*L*yz,   -14*L*yt,   3*L**2*yz,     -3*L**2*yy ],
651    ])
652
653    mat[6:12,6:12] = np.array([
654        [0,         0,          0,          0,          0,              0               ],
655        [0,         156*yy,    156*yz,    147*yt,    22*L*yz,      -22*L*yy        ],
656        [0,         156*zy,    156*zz,    147*zt,    22*L*zz,      -22*L*zy        ],
657        [0,         147*ty,    147*tz,    140*tt,    21*L*tz,      -21*L*ty        ],
658        [0,         22*L*zy,   22*L*zz,   21*L*zt,    4*L**2*zz,   -4*L**2*zy      ],
659        [0,         -22*L*yy,   -22*L*yz,   -21*L*yt,   -4*L**2*yz,    4*L**2*yy   ],
660    ])
661
662    return mat * L / 420
def bar_foundation_stiffness(L, kx, ky, kz):
665def bar_foundation_stiffness(L, kx, ky, kz):    #only z and y currently, will be extended!
666    mat = np.vstack([ 
667        [kx*1/4, 0, 0, 0, 0, 0,     kx*1/4, 0, 0, 0, 0, 0],     #ux1
668        [0, ky*1/3, 0, 0, 0, 0,     0, ky*1/6, 0, 0, 0, 0],     #uy1
669        [0, 0, kz*1/3, 0, 0, 0,     0, 0, kz*1/6, 0, 0, 0],     #uz1
670
671        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #rx1
672        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #ry1
673        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #rz1
674
675        [kx*1/4, 0, 0, 0, 0, 0,     kx*1/4, 0, 0, 0, 0, 0],     #ux2
676        [0, ky*1/6, 0, 0, 0, 0,     0, ky*1/3, 0, 0, 0, 0],     #uy2
677        [0, 0, kz*1/6, 0, 0, 0,     0, 0, kz*1/3, 0, 0, 0],     #uz2
678
679        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #rx2
680        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0],          #ry2
681        [0, 0, 0, 0, 0, 0,          0, 0, 0, 0, 0, 0]]) * L      #rz2
682    return mat