2. Static Analysis

2.1. Cable Truss

Location of example: $B2EXAMPLES/static/cable_truss

A truss composed of prestressed cables deforms under its own weight (simulated here by forces concentrated at the truss free nodes) and under an additional single force applied at one node. The model has been presented by Kar et.al [kar] . The non-linear analysis can be executed in one or in several stages being executed one after the other in the order they are enumerated in the case definition. Each stage references a case specifying sets of boundary conditions and analysis parameters and the solution will be incremented to a load or time parameter of 1.0 unless specified otherwise by the sfactor option of the stage.

_images/examples.static.cabletruss.mesh.png

Cable truss mesh. Nodes are marked branch:node, where branch is the branch identifier and node the external node identifier. Boundary conditions (red) are marked LLL, i.e displacements locked in all 3 directions at nodes 1, 2, 3, 6, 7, 10, 11, 12).

Two different cases are presented, all three leading to the same result. The problem must be solved in 3 steps (see also note below):

  1. Run the input processor b2ip++, defining all 3 cases (see input file cable.mdl.

  2. Run for case 1 with the solver option adir cases=1.

  3. Run for case 2 with the solver option adir cases=2.

The database now contains the 3 runs and the load-stress curves can be extracted (see diagrams below).

Note

The non-linear B2000++ solver can treat more than one nonlinear case in a single model. However, the cases must be executed one by one.

Case 1: The analysis is performed in one stage, with all boundary conditions (dead weight plus load) applied simultaneously. Case 1 tries to solve for all boundary conditions in one stage. Essential boundary conditions are defined in set 1 (option ebc 1). Loads are defined in natural boundary condition sets 1 and 2. The solver will try to solve for the final default load scale factor of 1.0 (when no stage is defined). Step control: In case the step size is too large, steps can become smaller, but not less that 0.01 (option step_size_min 0.01). :

case 1
    title 'Dead weight and useful loads in 1 set'
    analysis nonlinear
    ebc 1  # Zero constraints
    nbc 1  # Dead weight
    nbc 2  # Useful load
    step_size_init 0.1
    step_size_min 0.01
    step_size_max 0.1
    gradients 1
end

The output from B2000++ is listed below (it is obtained with the B2000+ command option -l):

b2000++ -l 'info of all in err' -adir 'cases=1' b2test.b2m
INFO:all:10:45:58.554: Start
INFO:solver.static_nonlinear:10:45:58.555: Start static nonlinear solver for case 1
INFO:domain:10:45:58.555: Total number of DOfs: 36.
INFO:solver.static_nonlinear:10:45:58.556: Start prediction of increment 1 of stage id 1 for time 0.
INFO:solver.static_nonlinear:10:45:58.591: Start correction of increment 1 of stage id 1 for time 1.
INFO:solver.static_nonlinear:10:45:58.595: End static nonlinear analysis for case 1 with success.
INFO:all:10:45:58.595: End of execution

Case 2: Case 2 solves the stages 21 and 21 in the order they are enumerated in the case 2 block. Again, essential boundary conditions are defined in ebc set 1 (option ebc 1). The stage with the dead weight, i.e case 21, refers to the natural boundary condition set 21 (option nbc 21), which is solved up to the default load scale factor 1.0. The solver then continues with the next stage, i.e. case 22, which refers to the natural boundary condition set 22 (option nbc 22), solving for the final load multiplier value of that stage, i.e. (default) 1.0. 22):

case 2
    title 'Stage 1 (own weight) and stage 2 (useful load)'
    analysis nonlinear
    stage 21 sfactor 1  # Dead weight
    stage 22 sfactor 1  # Useful load
    gradients 1
end

stage 21
    title 'Stage 1 (first stage, dead weight)'
    ebc 1  # Zero constraints
    nbc 1  # Dead weight
    step_size_init 0.1
    step_size_min 0.01
    step_size_max 0.1
end

stage 22
    title 'Stage 2 (useful load)'
    nbc 2  # Useful load added to all nbc's
    step_size_init 0.1
    step_size_min 0.01
    step_size_max 0.1
end

Both cases give the same result for the final load factor, although the path to the solution is not the same, i.e it depends on the load stages. Note that case 2 goes to a final load factor of 2 because of the 2 stages, while cases 1 and 2 go to the final load factor of 1.0.

_images/examples.static.cabletruss.case1.svg

Load case 1: Displacement and rod stress.

_images/examples.static.cabletruss.case2.svg

Load case 2: Displacement and rod stress. Note: The path is not the same as for case 1.

Deformations can be viewed with baspl++, the viewer script view.py producing the deformation plot with the cable stresses shown below for case 1. To run the viewer from the shell, type

baspl++ -t view.py
_images/examples.static.cabletruss.deformed.png

Deformed cable truss shape and cable stress (colored elements) and applied force.

In addition, node displacements of all cases can be printed with the B2000++ model browser b2browser.

b2browser cable

and click Solutions, Case 1, Displacements, Cycle 10 (the final cycle):

Displacements and Rotations (DISP), dataset "DISP.1.10.0.1", cycle=10, subcycle=0, case=1

      EID        UX           UY           UZ           RX           RY           RZ       Amplitude
        1 G           0            0            0            0
        2 G           0            0            0            0
        3 G           0            0            0            0
        4 G     -0.3509      -0.3509        0.887        1.016
        5 G     -0.5581        3.269       -2.816        4.351
        6 G           0            0            0            0
        7 G           0            0            0            0
        8 G       3.269      -0.5581       -2.816        4.351
        9 G       3.615        3.615        16.31        17.09
       10 G           0            0            0            0
       11 G           0            0            0            0
       12 G           0            0            0            0
Largest amplitude=17.0906
Absolute printout cutoff value threshold=1e-100

2.2. Cylinder with RBE/FMDE Elements

Location of example: $B2EXAMPLES/static/cylinder_rbe

This example demonstrates RBE Elements; RBE (rigid body element) example and FMDE Elements; FMDE (Force-Moment Distribution Element) example elements for defining boundary conditions (constraints) applied to a a cylinder section Boundary conditions; Apply with RBE elements. The model is a cylinder clamped at one end \(z=0\) and loaded at the other end \(z=L\). To introduce the corresponding boundary conditions, RBE/FMDE elements are introduced:

  • The clamping condition is formulated with a RBE/FMDE element connecting all DOFs of a node at (0,0,0) to all circumferential nodes of the cylinder at \(z=0\). The constraints (mdl.ebc) will be applied to that node. The RBE element acts as a rigid diaphragm in the xy-plane.

  • The load condition is formulated with a RBE element connecting a node at (0,0,L) with all circumferential nodes of the cylinder at \(z=L\). The force (nbc) will be applied to that node. The RBE/FMDE element acts as a rigid diaphragm in the xy-plane.

_images/examples.static.cylinderrbe.mesh.png

Cylinder FE mesh: Q9 shell element mesh (blue) and RBE/FMDE mesh (red).

Dimensions and material parameters of the cylinder model:

L = 1.0    # Length of cylinder [m]
R = 0.3    # Radius of cylinder [m]
T = 0.0005 # Thickness of cylinder [m]
Fx = 1.e5  # Axial load [N]
E = 72.4e9 # Young's modulus [N/m2]
_images/examples.static.cylinderrbe.disp.png

Cylinder model: Deformation amplitude, traction case.

_images/examples.static.cylinderrbe.disp-z.png

Cylinder model: Deformation amplitude of RBE/FMDE element, traction case.

2.3. Laminate Material (1 Element Test)

Location of example: $B2EXAMPLES/static/laminate1

This one-element model consists of a rectangular four-node Q4 shell element in the global x-y plane The model is “statically determinate”, i.e the number of constrained degrees of freedom is equal to the number of rigid body movements, i.e 3 translations and 3 rotations. The plate is pulled in the element x-direction (fiber 1 direction, see figure below).

The purpose of this example is to

  • Demonstrate laminate material in shells with a one layer material with 0°, 90° fiber orientations, and with a four layer 0°/90°/90°/0° material.

  • Explain the basic concepts of coordinate systems.

Forces (natural boundary conditions) are added as point forces in order to illustrate the concept of node-local coordinate systems. The total force applied is selected such that the stress in fiber direction 1 becomes 500 MPa:

\[F_{tot}=\sigma_{11}*B*T = 500*2.3*0.15 = 194 N\]

The first test 0.mdl defines a model in the branch-global (=global-global) x-y plane:

  • The branch-global system is identical to the global-global system (no branch orientation specified).

  • The coordinates are formulated in the branch-global system.

  • The element local x-axis coincides with the branch-global x-axis.

  • The node-local directions of the degrees-of-freedom coincide with the branch-global direction.

  • The zero constraints (essential boundary conditions, ebc) do not need node-local transformations.

  • The applied forces \(F=0.5*F_{tot}\) in the branch-global x-direction do not need node-local transformations.

_images/examples.static.laminate1.mesh.svg

One-element laminate test: Element local x-axis coincides with branch-global x-axis. Gray: mesh, blue: deformed mesh, red: constrains, green: forces.

The second test 30.mdl defines a model in the branch-global (=global-global) x-y plane with the element local x-axis rotated by 30° around the branch-global z-axis with respect to the branch-global x-axis:

  • The branch-global system is identical to the global-global system (no branch orientation specified).

  • The coordinates are formulated in the branch-global system.

  • The orientation of the values at the mesh nodes (forces, constraint) are formulated in a node-local orientation system with the MDL transformation command

    transformations
        1 cartesian 0 0 0  0 0 1  (cos30) (sin30) 0
    end
    

    and referred to by the dof_ref parameter of the MDL nodes command:

    nodes
        dof_ref 1 # refers to transformation 1
        1 0.0 0.0 0.0
        2 7.53442101 4.35 0.0
        3 6.38442101 6.34185843 0.0
        4 -1.15 1.99185843 0.0
    end
    
  • The zero constraints (essential boundary conditions, ebc) are formulated in the node-local system, because the nodes involved have a local transformation.

  • The applied forces \(F=0.5*F_{tot}\) are formulated in the node-local system, because the nodes involved have a local transformation.

_images/examples.static.laminate1.mesh30.svg

One-element laminate test: Model rotated 30° around the branch-global z-axis. Gray: mesh, blue: deformed mesh, red: constrains, green: forces.

The third test 30btrf.mdl defines a model in a branch-global

transformed with respect to the global-global system by a branch_orientation command:

  • The branch-global system is rotated with the MDL branch_orientation command by 30° around the global-global z-axis:

    branch_orientation
        rotate axis z angle 30
    end
    
  • The coordinates are formulated in the branch-global system and are identical to the coordinates of model 0.mdl.

  • The element local x-axis coincides with the branch-global x-axis.

  • The node-local directions of the degrees-of-freedom coincide with the branch-global direction (no node-local transformation).

  • The zero constraints (essential boundary conditions, ebc) do not need node-local transformations.

  • The applied forces \(F=0.5*F_{tot}\) in the branch-global x-direction do not need node-local transformations.

Comparison of the results of the three models (the solutions of the analysis) is discussed. To resume, please note that:

  • DOF fields, i.e. displacements and rotations, are stored on the B2000++ model database with respect to the branch-global coordinate system. Values at nodes with node-local coordinate systems (orientations) are transformed to the branch-global orientation.

  • Element stress and strain tensors as well as any other sampling point field values, i.e element-related derived quantities, are stored on the B2000++ model database with respect to the branch-global coordinate system.

Displacements of nodes 2 and 3 as stored on the B2000++ model database are compared in the table below:

Comparison of displacements
   Model  EID           DX           DY    Amplitude
     "0"    2      0.02962            0      0.02962
            3      0.02962    -0.002349      0.02971

    "30"    2      0.02565      0.01481      0.02962
            3      0.02683      0.01278      0.02971

"30btrf"    2      0.02962            0      0.02962
            3      0.02962    -0.002349      0.02971

Displacements of model 0 and 30btrf are identical, because both models are formulated in the same branch global configuration.

Stresses at specific points can be plotted or printed with the b2plot_stress_stack script.

script. Example: Plot stresses through the shell thickness z' at integration point 1 of element 1.

b2plot_stress_stack -e 1 -p 2 --case 1 --print -f STRESS 30;
_images/examples.static.laminate1.STRESS-1-2.svg

Model 30: Stresses through shell thickness z' at integration point 1 of element 1.

Print stresses through the shell thickness z' at integration point 1 of element 1.

b2plot_stress_stack -e 1 -p 2 --case 1 -f STRESS --print 30

Layer          z   Sigma_xx   Sigma_yy   Sigma_zz   Sigma_xy   Sigma_yz   Sigma_xz
    0      -0.08        701        246   8.98e-35        394   3.56e-17  -2.06e-17
    0      -0.04        701        246   2.25e-35        394   1.78e-17  -1.03e-17
    1      -0.02       48.7        3.9  -3.01e-36       38.8  -6.57e-18   3.79e-18
    2       0.02       48.7        3.9  -3.01e-36       38.8   6.57e-18  -3.79e-18
    3       0.04        701        246   2.25e-35        394  -1.78e-17   1.03e-17
    3       0.08        701        246   8.98e-35        394  -3.56e-17   2.06e-17

2.4. Mixed Element Types Stress Tests

Location of example: $B2EXAMPLES/static/mixed

Checks and demonstrates the availability of stresses in a FE model with beam, rod and shell elements. The B2000++ solver computes and stores

  • Stresses in shells and solid elements.

  • Forces and moments in beam elements. Stresses in beam elements must be computed in a post-processing step with Simples.

  • Forces and stresses/strains in rod elements.

The model consist of a combination of shell and beam elements forming a steel plate on 4 steel columns and diagonals (see FE model plot below) and it is generated with the Simples modeler (see mixed.py). mixed.py also launches the solver.

_images/examples.static.mixed.mesh.png

Mixed problem: Mesh generated with the Simples modeler and viewed with the baspl++ script bv-view-mesh.py).

The plate is loaded with a pressure load of 5000 Pa. The columns are clamped at the base.

The MDL input file mixed.mdl is generated with the Simples model generation script mixed.py. mixed.py also launches the solver.:

python3 mixed.py
mixed.py also launches the solver, producings the model

database mixed.b2m.

The post-processor (baspl++) scripts bv-view-mesh.py display the mesh and bv-view-deformed.py the deformed shape:

_images/examples.static.mixed.deformed.png

Mixed problem: Deformed shape.

In the current version of B2000++ (4.4) beam stresses cannot be computed with the B2000++ solver. This is illustrated in the figure below displaying the failure index: Beam elements are gray (no failure index found).

_images/examples.static.mixed.findex.png

Mixed problem: Failure index.

Beam stresses must be evaluated in a post-processing step with Simples. The stress extraction script print_beam_rod_stresses.py illustrates this step: It computes beam section stresses and prints them for one of the columns of the model. It also prints the maximum of \(\sigma_{xx}\) (min and other components are extracted similarly with the get_max() and get_min() functions).

"""Simples script to extract and print beam and rod stresses from
"mixed" model."""

import sys
import numpy
import simples as si

def print_beam_stresses(model, elset=None, branch=1, case=1, cycle=0):
    """Computes and prints beam section stresses for all elements of list
    ``elset``.

    """
    field = si.BeamStressField(model, branch=branch, case=case, cycle=cycle)
    if field is None:
        print(f"Beam stress field not found")
        return
    for r in field.get_stresses(elset):
        print(f"Beam stresses, element {r.eid}")
        print(f"   {'ip':>2} {'i':>2} {'sxx':^10} {'st':^10}")
        for ip, stress in enumerate(r.stresses):
            for i, s in enumerate(stress):
                print(f"   {ip:2d} {i:2d} {s[0]:10.4g} {s[1]:10.4g}")

def print_rod_stresses(model, elset=None, branch=1, case=1, cycle=0):
    """Computes and prints rod axial stresses for all elements of list
    ``elset``.

    """
    field = si.SamplingPointField(model, 'STRESS_SECTION_ROD', branch=branch, case=case, cycle=cycle)
    if field is None:
        print(f"Rod stress field not found")
        return
    if elset is not None:
        _eids = tolist(elset)
    else:
        _eids = model.get_branch().get_elem_eids()
    for eid in _eids:
        elstresses = field.get(eid)
        if elstresses is not None:
            print(f"Rod stresses, element {eid}\n   {'ip':>2} {'sxx':^10}")
            for ip, stress in enumerate(elstresses):
                 print(f"   {ip:2d} {stress[0]:10.4g}")

dbname = ('mixed.b2m' if len(sys.argv) < 2 else sys.argv[1])

model = si.Model(dbname)

# Compute beam stresses for the first column "column1". Columns 2-4
# idem.

elset = model.get_elem_set("column1")
print_beam_stresses(model) # prints all element stresses
#-demo :print_beam_stresses(model, elset) # prints all beam element stresses of elset
print_rod_stresses(model) # prints all rod element stresses of elset

Executing the script with Python3:

python3 print_stresses.py

Beam stresses, element 1
   ip  i    sxx         st
    0  0 -1.665e+07          0
    0  1 -1.131e+07          0
    0  2 -1.665e+07          0
    0  3 -1.382e+07          0
Beam stresses, element 2
   ip  i    sxx         st
    0  0 -1.195e+07          0
    0  1 -1.275e+07          0
    0  2 -1.195e+07          0
    0  3 -1.238e+07          0
Beam stresses, element 3
   ip  i    sxx         st
    0  0  -7.24e+06          0
    0  1  -1.42e+07          0
    0  2  -7.24e+06          0
    0  3 -1.093e+07          0
Beam stresses, element 4
   ip  i    sxx         st
    0  0 -2.536e+06          0
    0  1 -1.564e+07          0
    0  2 -2.536e+06          0
    0  3 -9.487e+06          0
Stress max(Sxx) in element smax[0], sxx= 1.051e+06

2.5. Node-Local Transformation

Location of example: B2EXAMPLES/static/node_local_transformation

This basic test explains the local orientations of values at nodes, such as boundary condition or DOF values. The test model is a cantilever beam modeled with 4 B2.S.RS beam elements. It is clamped at node 1 and loaded with a force in the global y direction at node 5. Three nodes are defined in a local system by defining transformations

transformations
  # Transformation 3: Rotation around z-axis by -90 deg
  3 cartesian  0. 0. 0.  0. 0. 1.  0. -1 0. # x_local = -y

  # Transformation 5: Rotation around z-axis by 90 deg, x_local = y: Ux_local = Uy_global
  5 cartesian  0. 0. 0.  0. 0. 1.  0. 1. 0.
end

and assigning them to nodes

nodes
  1 0. 0. 0.
  2 2.5 0. 0.
  transformation 3
  3 5. 0. 0.
  transformation 0 # Restore global
  4 7.5 0. 0.
  transformation 5
  5 10. 0. 0.
end

The node-local Cartesian coordinate systems of boundary conditions and DOF values of the mesh nodes are as follows:

  • Node 1 has (default) global orientation (no transformation).

  • Node 2 DOF values are defined in the local system of transformation 3. Since no boundary conditions are define at node 3 no transformation will be applied by the solver.

  • Node 4 has (default) global orientation (no transformation).

  • Node 5 has the assigned local transformation 5: DOF 1 points to the global direction y, DOF 2 to global direction -x, and DOF 3 to global direction z. Transformations will be applied by the solver.

_images/examples.static.node_local_transformation.node_local.png

Beam mesh (black). Nodes marked by black squares. Local orthogonal coordinate systems (red: x-axis, green: y-axis, blue: z-axis). baspl++ script bv-node-local-coordinate-systems.py generates this plot.

nbc and ebc conditions must be formulated and input in the node-local system. Thus, defining an applied force (value=30) in the global y-direction with the nbc command

# Applied force Fy in global y-direction (DOF 2). Local: Fx (DOF 1)
nbc 1
  value 30 dof Fx nodes 5
end

requires to assign the force to DOF 1 due to the local transformation.

Results are by default generated and stored with respect to the (branch) global system. For node 5 the generated applied forces are stored on the database in the (branch) global system:

Applied forces at node 5

Node

Sys

Fx

Fy

Fz

Mx

My

Mz

5

G

0

30

0

0

0

0

Likewise, the resulting displacements and rotations are stored on the database in the (branch) global system:

Displacements and rotations at node 5

Node

Sys

Ux

Uy

Uz

Rx

Ry

Rz

5

G

0

1.078

0

0

0

0.15

2.6. Scordelis-Lo Roof (Cylindrical Shell)

Location of example: file:$B2EXAMPLES/static/scordelis_lo_roof

The Scordelis-Lo roof is a classical test case for linear static analysis named after the authors who first reported it [scordelis]. The test has been proposed as a standard FE test by MacNeal and Harder [macneal1985] and it can be found in the NAFEMS test case description.

The structure consist of a concrete cylindrical shell roof supported by rigid diaphragms on both curved sides (see figure below), and it is loaded by a surface traction of 90psf in the negative y-direction. The total length of the cylinder is 50ft, the radius of the cylinder 25ft, the half opening angle 40°, and the thickness 0.25 ft. The material properties are isotropic, with the elastic modulus of 4.32e8psf and the Poisson number set to 0 or, in some references, to 0.17. Von Mises failure stress is 300000psf.

_images/examples.static.scordelis.mesh.png

Scordelis-Lo Roof: Mesh.

The diaphragms are defined such that a simply supported boundary condition results: Displacements in x- and y-directions are suppressed, the suppression of displacements in the z-direction being assumed by the symmetry condition along \(P_4\)-\(P_1\).

The test compares the linear static solution, specifically the displacement in the y-direction of point \(P_4\), against known solutions (the reference solution for point \(P_4\) is -0.3024).

The solution of this problem with B2000++ consist of the following steps:

  • Create the finite element model with the MDL language describing the problem in a text file.

  • Solve the problem by submitting the text file to the B2000++ solver. The solver creates a B2000++ database.

  • Examine the results by reading the B2000++ database.

A finite element model for B2000++ is created by describing the problem geometry, creating a finite element mesh, adding boundary conditions, and adding solver instructions. The complete model can be formulated with the B2000++ Model Description Language (MDL). Alternatively, parts of the model, such as the mesh and the boundary conditions, can be produced with an external program and brought in the B2000++ input file.

The Scordelis-Lo roof problem can easily be produced with the B2000++ Model Description Language. The necessary steps will now be explained. All steps result in instructions added with an editor to a single B2000++ input file.

Geometry and mesh: The cylinder geometry and the mesh produced by the geometry can be described with an epatch (element patch). Element patches allow for generating simple building blocks, such as plates, cylinders, etc. The Scordelis-Lo cylindrical roof is described by an epatch of type cylinder, which implies the creation of a surface and a surface mesh:

epatch 1
    geometry cylinder
    type (ELTYPE)     # Element typ
    ne1 (NE) ne2 (NE) # N. of elements along sides
    thickness 0.25    # Thickness of shell
    mid 1             # Material identifier
    phi1 50.0         # Opening angle of half cylinder
    phi2 90.0         # Closing angle of half cylinder
    radius 25.0       # Radius of cylinder
    length 25         # Half length of cylinder
    # local edges means that all DOF's of the nodes along the edges are
    # expressed in local, i.e. cylindrical, coordinates.
    if (local == 1) {local edges}
end

Note that the epatch must contain the element type (eltype keyword), the number of elements to be generated along each side of the patch (ne1 and ne2 keywords), and material reference (mid keyword). The definition of the epatch is further illustrated in the figure below, with the position of the vertices P1 to P4 and the edges E1 to E4.

In this particular configuration the local and global orientations are the same, i.e. no need to define cylindrical coordinates, unless the resulting displacements, rotations, forces, and moments should be expressed in cylindrical coordinates (see instruction if (local == 1) {local edges} above).

_images/examples.static.scordelis.epatch.png

Scordelis-Lo Roof: epatch vertex points and edge definitions. See also the generated mesh.

Material: The material is the same as the one in the original problem description. The definition is straightforward with the material'' command, which identifies the material with id ``1

material 1
    type isotropic
    e 4.32e8       # Young's modulus
    nu 0.          # Poisson number
    density 90
    failure von_mises
      r 300000
    end
end

Boundary conditions: Since all boundary conditions apply either to edges of the epatch or to the epatch surface they can be generated for all nodes of an edge by referring to the epatch edges which in their turn are automatically generated. Essential boundary conditions (here: displacement constraints) are defined with the ebc command:

ebc 1
    value 0 # Values to be assigned to DOF's
    # Diaphragm = Supported edge (global=local orientation) along
    # edge e1 from P1 to P2
    dof [UX UY] epatch 1 e1
    # Symmetry condition along edge e2 from P2 to P3
    dof [UX RY RZ] epatch 1 e2
    # Symmetry condition along edge e3 from P3 to P4
    dof [UZ RX RY] epatch 1 e3
end

The figure below shows the essential boundary conditions as defined above: At each node of the edges a pattern XXXXXXX is plotted, where the position of X in the pattern designates the degree of freedom (DOF) to be affected and X taking the values F (free, i.e. unconstrained), L (locked, i.e. constrained to 0), or V (constrained to a value).

_images/examples.static.scordelis.bc.png

Scordelis-Lo Roof: Essential boundary conditions generated by ebc 1.

The surface pressure loads are applied with the natural boundary condition command nbc. Note that the reference frame for the surface tractions is the branch reference frame and thus has to be specified explicitly with the system option, because the default reference frame is the local reference frame relative to the element surface.

nbc 1 type surface_tractions
    system branch # Branch referential (default is local!)
    surface_tractions 0 -90 0
    epatch 1 f7   # Apply to mid-surface of shell!
end

What remains now to be specified are (1) definition of the case, comprising the type of analysis (2) the boundary condition sets to be applied to the case (3) the analysis directives, i.e. the cases to be solved for:

case 1
    analysis linear
    nbc 1
    ebc 1
end

adir
    case 1
end

The analysis is run with the shell command:

b2000++ scordelis_lo_roof.mdl

which generates the B2000++ database scordelis_lo_roof.b2m. Specifying the log option to the b2000++ command prints a short summary of the steps which the solver performs:

$ b2000++ -l info scordelis_lo_roof.mdl
INFO:all:13:57:24.656: Start
INFO:solver.static_linear:13:57:24.658: Start static linear solver
INFO:domain:13:57:24.658: Total number of DOfs: 422.
INFO:solver.static_linear:13:57:24.659: Element matrix assembly
INFO:solver.static_linear:13:57:24.660: Resolve the linear problem
INFO:solver.static_linear:13:57:24.664: Compute gradients and reaction forces
INFO:solver.static_linear:13:57:24.664: End of static linear solver
INFO:all:13:57:24.664: End of execution

Results

The summary of the analysis results for the 8 by 8 element mesh is displayed in the table below for element type Q4.S.MITC.E4:

Q4.S.MITC.E4 results (8 by 8 element mesh)

What

Reported

B2000++

Error (%)

Displacement DY at point P3 [ft]

-0.3024

-0.304196

0.6

Results can be examined with several B2000++ tools, such as the baspl++ post-processor, the Simples scripting tools, the *b2print_ text output tools, or the B2000++ data browser. We start with the baspl++ viewer to get an image of deformed shape of the structure. The figure below is obtained by launching baspl++ from the shell with the viewer script view.py:

baspl++ -t plot_displacements.py

and shows the deformed structure (colored) and the initial undeformed one (black dashed grid). Please read the file view.py.

_images/examples.static.scordelis.disp.png

Scordelis-Lo Roof: Displacement field.

baspl++ -t plot_stresses.py

plots the per element maximum von Mises comparison stress:

_images/examples.static.scordelis.vmises1.png

Scordelis-Lo Roof: von Mises stresses (max of element).

baspl++ -t plot_stresses_sampling.py

plots the von Mises comparison stress at each integration (sampling) point in the form of spheres, where color and sphere size are proportional to the intensity:

_images/examples.static.scordelis.vmises2.png

Scordelis-Lo Roof: von Mises stresses (all sampling points).

Values can be extracted with Simples scripts. The script print_displacement.py extracts the displacement at point P3 and prints it:

python print_displacement.py

producing the output

Displacement at point P3 (node 273): DY=-0.301733
Reference solution DYREF=-0.3024
Error=-0.22%

and the script print_stresses.py extracts and prints stress values at selected elements prints it:

python3 print_vmises.py

Max failure index fmax=1.124541, element=241

Data can be viewed with the B2000++ data browser. Example: View the failure index:

_images/examples.static.scordelis.b2browser.png

Scordelis-Lo roof: View failure index values with the B2000++ browser.

A convergence study increasing the number of cells (mesh density) with Simples scripting is very easily performed (script convergence.py). The graph produced by the script is displayed below:

_images/examples.static.scordelis.convergence_all.svg

Scordelis-Lo roof: Convergence graph for selected element types.

The baspl++ viewer allows for off-screen (batch mode) rendering. Batch mode is particularly useful if images are to be generated on a plaform without graphic rendering. Batch mode must be launched with a script. Example script plot_displacements_batch.sh:

baspl++ --batch --mesa << /

#!/bin/sh
# Example script to generate a plot "output.png" with baspl++ in batch
# (off-screen rendering) mode.
# Observe s.add(p) below: This is done automatically in non batch
# mode, but required in batch scrpting mode. Failing to add the part
# to the scene will produce an empty (scene) plot!

s=Scene()
s.size = (1200, 900)
s.modelview.matrix = [
    [ 0.03505605,  0.00484255,  0.01534965,  0.00000000],
    [-0.00926785,  0.03614971,  0.00976166,  0.00000000],
    [-0.01315936, -0.01255918,  0.03401601,  0.00000000],
    [ 0.08682165, -0.65933049, -1.25635135,  1.00000000],
]
s.scene_axes_show = False
m = Model('scordelis_lo_roof.b2m')
p = NPart(m)
s.add(p)  # Add part to scene
p.deform.scale = [5.0, 5.0, 5.0]
p.edge.show = True
p.elements.extract = True
f = Field(m, 'DISP', case=1)
p.deform.field = f
p.contour.compname = '2'
p.contour.colmap_default.texture.mode = 'step'
p.contour.field = f
p2 = NPart(m)
p2.edge.show = True
p2.edge.width = 3.0
p2.elements.extract = True
p2.face.show = False
s.export_to_file('output.png')
/
_images/examples.static.scordelis.off-screen-rendering.png

File “output.png” produced in batch mode by script “plot_displacements_batch.sh”.

2.7. Linear Analysis with Subcases

Location of example: $B2EXAMPLES/static/subcases

Th example illustrates the concept of subcases (see case definition). During linear analysis, subcases can be applied to a common set of essential boundary conditions, thus reducing solution time, because the constraints are applied once and the global matrix is factored once. Each subcase is a new right-hand side and needs back-substitution only.

The model is the same as the Buckling of Composite Plate example: A simply-supported thin quadratic laminated plate of 500 mm by 500 mm is loaded with 2 load cases and same essential boundary conditions.

The first load case (case 100) constrains the plate to simply supported condition (essential boundary condition 1). Under this condition, 2 natural boundary conditions (nbc) force sets are applied as subcases:

case 100
    analysis linear
    gradients 1
    ebc 1
    subcase 1 # Subcase 1: Pressure load (nbc 1)
        nbc 1
    end
    subcase 2 # Point load applied at the plate's mid-point (nbc 2)
        nbc 2
    end
end

The second load case (case 200) is the same as case 100 and serves as test for subcase reference.

case 100
    analysis linear
    gradients 1
    ebc 1
    subcase 1 # Pressure load (nbc 1)
        nbc 1
    end
    subcase 2 # Point load at mid-node
        nbc 2
    end
end
_images/examples.static.subcases.definition.png

Geometry and essential boundary conditions of plate.

Note that B2000++ always solves for all subcases of a case, i.e specifying

adir
    cases [100 200]
end

will solve for subcase 1 and subcase 2 of cases 100 and 200.

The input file plate.mdl contains all definitions for solving the 2 cases with 2 subcases each. The essential and natural boundary condition definitions are not listed here (see MDL input file plate.mdl).

The deformed shapes of the subcases are displayed in the figures below. Note that the selection of the subcase in baspl++ and Simples is obtained with the subcase= or mode= parameter.

_images/examples.static.subcases.disp-1-1.png

Linear Analysis of Plate with subcases: Amplified deformed shape, case 1, subcase 1.

_images/examples.static.subcases.disp-1-2.png

Linear Analysis of Plate with subcases: Amplified deformed shape, case 1, subcase 2.

_images/examples.static.subcases.disp-2-1.png

Linear Analysis of Plate with subcases: Amplified deformed shape, case 2, subcase 1.

_images/examples.static.subcases.disp-2-2.png

Linear Analysis of Plate with subcases: Amplified deformed shape, case 2, subcase 2.

Execute the test with

b2000++ plate.mdl

Solutions can the be examines with b2print_solutions``. Example: Print the displacements and rotations of mid-node 128 for case 200 and subcase 2:

b2print_solutions --case=200 --subcase=2 --item=128  plate.b2m/

DOF Field "DISP", mesh=1, case=200, cycle=0, subcase=2, selection="Selected nodes"
EID           Ux           Uy           Uz           Rx           Ry           Rz
128     2.24e-21    -1.52e-20        2.128     0.001757     -0.01189            0

Plot the deformed shape with baspl++ script bv-deform.py.

Finally, the test runner executes the test and heck the force and displacement results (see __init__.py. Execute from the test directory with the shell:

b3tesrunner .