p_shell#

Purpose

Element property function for shells and plates (flat shells)

Syntax

il = p_shell('default');
il = p_shell('database ProId name'); 
il = p_shell('dbval ProId name');
il = p_shell('dbval -unit TM ProId name');
il = p_shell('dbval -punit TM ProId name');
il = p_shell('SetDrill 0',il);

Description

This help starts by describing the main commands : p_shell Database and Dbval. Supported p_shell subtypes and their formats are then described.

Database,Dbval, ...#

p_shell contains a number of defaults obtained with the database and dbval commands which respectively return a structure or an element property row. You can select a particular entry of the database with using a name matching the database entries.

You can also automatically compute the properties of standard shells with

kirchhoff eKirchhoff shell of thickness e (is not implemented for formulation 5, see each element for available choices)
mindlin eMindlin shell of thickness e (see each element for choices).
laminate MatIdi Ti ThetaiSpecification of a laminate property by giving the different ply MatId, thickness and angle. By default the z values are counted from -thick/2, you can specify another value with a z0.

You can append a string option of the form -f i to select the appropriate shell formulation. The different formulations are described under each element topology (tria3,   4, ...) For example, you will obtain the element property row with ProId 100 associated with a .1 thick Kirchhoff shell (with formulation 5) or the corresponding Mindlin plate use

 il = p_shell('database 100 MindLin .1')
 il = p_shell('dbval 100 kirchhoff .1 -f5')
 il = p_shell('dbval 100 laminate z0=-2e-3 110 3e-3 30 110 3e-3 -30')
 il = fe_mat('convert SITM',il);
 il = p_shell(il,'dbval -unit TM 2 MindLin .1') % set in TM, provide data in SI
 il = p_shell(il,'dbval -punit TM 2 MindLin 100') % set in TM, provide data in TM

For laminates, you specify for each ply the MatId, thickness and angle.

Shell format description and subtypes#

Element properties are described by the row of an element property matrix or a data structure with an .il field containing this row (see section 7.4). Element property functions such as p_shell support graphical editing of properties and a database of standard properties.

For a tutorial on material/element property handling see section 4.5.1. For a reference on formats used to describe element properties see section 7.4.

p_shell currently only supports two subtypes

1 : standard isotropic#

  [ProID type   f d O   h   k   MID2 RatI12_T3 MID3 NSM Z1 Z2 MID4]
type identifier obtained with fe_mat('p_shell','SI',1).
f 0 use default of element. For other formulations the specific help for each element (quad4, tria3, ...), each formulation specifies integration rule.
d-1no drilling stiffness. The element DOFs are the standard translations and rotations at all nodes (DOFs .01 to .06). The drill DOF (rotation .06 for a plate in the xy plane) has no stiffness and is thus eliminated by fe_mk if it corresponds to a global DOF direction. The default is d=1 (d is set to 1 for a declared value of zero).
 darbitrary drilling stiffness with value proportional to d is added. This stiffness is often needed in shell problems but may lead to numerical conditioning problems if the stiffness value is very different from other physical stiffness values. Start with a value of 1. Use il=p_shell('SetDrill d',il) to set to d the drilling stiffness of all p_shell subtype 1 rows of the property matrix il.
h plate thickness.
kkshear correction factor (default 5/6, default used if k is zero). This correction is not used for formulations based on triangles since tria3 is a thin plate element.
RatI12_T3 Ratio of bending moment of inertia to nominal T3/I12 (default 1).
NSM Non structural mass per unit area.
MID2 material property for bending. Defauts to element MatId if equal to 0.
MID3 material property for transverse shear.
z1,z2 (unused) offset for fiber computations.
MID4 material property for membrane/bending coupling.

Shell strain is defined by the membrane, curvature and transverse shear (display with p_shell('ConstShell')).

 
    (9.25)

2 : composite#

  [ProID type   Z0 NSM SB FT TREF GE LAM MatId1 T1 Theta1 SOUT1 ...]
ProIDSection property identification number.
typeIdentifier obtained with fe_mat('p_shell','SI',2).
Z0Distance from reference plate to bottom surface.
NSMNon structural mass per unit area.
SBAllowable shear stress of the bonding material.
FTFailure theory.
TREFReference temperature.
EtaHysteretic loss factor.
LAMLaminate type.
MatIdiMatId for ply i, see m_elastic 1, m_elastic 5, ...
TiThickness of ply i.
ThetaiOrientation of ply i.
SOUTiStress output request for ply i.

Note that this subtype is based on the format used by NASTRAN for PCOMP and the formulation used for each topology is discussed in each element (see quad4, tria3). You can use the DbvalLaminate commands to generate standard entries.

 
    (9.26)

setTheta#

When dealing with laminated plates, the classical approach uses a material orientation constant per element. OpenFEM also supports more advanced strategies with orientation defined at nodes but this is still poorly documented.

The material orientation is the reference for plies. Any angle defined in a laminate command is an additional rotation. In the example below, the element orientation is rotated 30 degrees, and the ply another 30. The fibers are thus oriented 60 degrees in the xy plane. Stresses are however given in the material orientation thus with a 30 degree rotation. Per ply output is not currently implemented.

The element-wise material angle is stored for each element. In column 7 for tria3, 8 for quad4, ... The setTheta command is a utility to ease the setting of these angles. By default, the orientation is done at element center. To use the mean orientation at nodes use command option -strategy 2.

model=ofdemos('composite');
model.il = p_shell('dbval 110 laminate 100 1 30'); % single ply

% Define material angle based on direction at element MAP=feutil('getnormalElt MAP -dir1',model); bas=basis('rotate',[],'rz=30;',1); MAP.normal=MAP.normal*reshape(bas(7:15),3,3)'; model=p_shell('setTheta',model,MAP);

% Obtain a MAP of material orientations MAP=feutil('getnormalElt MAP -dir1',model); feplot(model);fecom('showmap',MAP)

% Set elementwise material angles using directions given at nodes. % Here a global direction MAP=struct('normal',ones(size(model.Node,1),1)*bas(7:9), ... 'ID',model.Node(:,1),'opt',2); model=p_shell('setTheta',model,MAP);

% Using an analytic expression to define components of % material orientation vector at nodes data=struct('sel','groupall','dir',{{'x-0','y+.01',0}},'DOF',[.01;.02;.03]); model=p_shell('setTheta',model,data); MAP=feutil('getnormalElt MAP -dir1',model); feplot(model);fecom('showmap',MAP)

model=p_shell('setTheta',model,0) is used to reset the material orientation to zero.

Technically, shells use the of_mk('BuildNDN') rule 23 which generates a basis at each integration point. The first vector v1x,v1y,v1z is built in the direction of r lines and v2x,v2y,v2z is tangent to the surface and orthogonal to v1. When a InfoAtNode map provides v1x,v1y,v1z, this vector is projected (NEED TO VERIFY) onto the surface and v2 taken to be orthogonal.

See also

Section 4.5.1, section 7.4, fe_mat