This input guide is a (nonofficial) companion of MATLAB program AEPH from Anisotropic elasticity with Matlab, Chyanbin Hwu, Springer, 2021, where additional information can be found in § 3.1.4, §15.1.7 (for BEM) and Appendix E.

All the labels used in AEPH are numbers, i.e. all data in the input files must be numeric.

Parameters enclosed by brackets [ ] are optional.

• input_control.txt (global parameter)
• input_elastic.txt (elastic property)
• input_thermal.txt (thermal property)
• input_loadstr.txt (problem parameter)
  • Ltype format for boundary element methods
  • Ltype=1
  • Ltype=4
  • Ltype=7
  • Ltype=8
  • Ltype=411
  • Ltype=611
  • Ltype=614
  • Ltype=622
• input_variable.txt (output parameter)
• input_xn.txt (node, BEM)
• input_node1.txt (element, BEM)
• input_bc.txt (boundary condition, BEM)
• input_temp.txt (thermal boundary condition, BEM)

Nmat Dtype Ltype Otype Icheck E0 h0 eps
Etype Ttype Ptype Vtype
[Etype Ttype Ptype Vtype]
[...]

Those in the 1st line are global control parameters. The material property parameters start from the 2nd line. Each line is a different material.

Parameters enclosed by brackets [ ] are optional.

E.g.

1 2 611 2 0 1e9 1e-3 1e-6
1 0 0 0

Total number of materials.

Valid problem dimension labels are:

• 1, 11, 111, 112, 12, 121, 122 Generalized plane strain.
• 2, 21, 211, 212, 22, 221, 222 Generalized plane stress.
• 3, 31, 311, 312, 32, 321, 322 Coupled stretching-bending.
• 4, 14, 104 Three-dimensional.

Solution label. 0~99 are boundary element method, and Ltype > 100 are analytical solutions (section number in Hwu (2021)). A complete list for BEM Ltypes can be found in §15.1.7.

Output scope label. Used with input_variable.txt, see that section for more information. Valid labels are:

• 1, 11, 12, 13: Curve.
  Additional output for analytical solutions (Ltype > 100):
  • 11: $\sigma_{ss}$, hoop stress (hole)
  • 12: $\sigma_{ss}, \sigma_{nn}, \sigma_{sn}$, internal surface traction (inclusion/interface?)
  • 13: $\sigma_{nn}$, normal stress (contact?)
• 2: Area.
• 3: Discrete points.

Enable/disable internal check. If and only if Icheck=1, values of A, B, mu, N1, N2, N3, N are verified through alternative approach and some identities. In some cases, also verifies the solution. No matter the checks are passed or not, the results will be displayed in the command window.

Reference Young's modulus for the nondimensionalization of N2, N3, and N. Suggested value=1.0e9.

Reference thickness used for the nondimensionalization of Aij,Bij and Dij. Suggested value=1.0e-3.

Perturbation ratio and error tolerance, suggested value=1.0e-6.

The suggested values of E0, h0, and eps are taken from comments in Main.m.

Used with input_elastic#.txt, see that section for more informations. Valid elastic property labels:

• 0 No elastic properties from input_elastic#.txt. Often used for piezoelectric or magneto-electric-elastic materials.

• 1 Isotropic.

• 2 Orthotropic.

• 3 Anisotropic, values in input_elastic#.txt are Cij.

• 4 Anisotropic, values in input_elastic#.txt are Sij.

• 5 Unidirectional fiber-reinforced composite.

• 6 Composite laminate.

Used with input_thermal#.txt, see that section for more information. Valid thermal property labels:

• 0 No thermal properties.
• 1 Isotropic.
• 2 Orthotropic.
• 3 Anisotropic, values in input_thermal#.txt are $k_{ij}$ and $\beta_{ij}$.
• 4 Anisotropic, values in input_thermal#.txt are $k_{ij}$ and $\alpha_{ij}$.
• 5 Unidirectional fiber-reinforced composites.
• 6 Composite laminates.

Piezoelectric properties input label. Used with input_piezo#.txt, see Ch. 11 of Hwu (2021) for more information. Valid labels are:

• 0 No piezoelectric properties.
• 1, 2, 3, 4
• 5 Electro-elastic laminates.
• 11, 12, 13, 14, 15, 16, 17, 18, 19 Magneto-electro-elastic materials.

Viscoelastic properties input label. Valid labels are:

• 0 No viscoelastic properties.
• 1, 2 Isotropic.
• 3 Standard linear viscoelastic solids.
• 4 Prony series.

Here goes the elastic properties.

Each material has an input file; e.g. four materials require four input_elastic#.txt files. The number 1 in the file name is its material number (order of appearance in input_control.txt).

Its format depends on Etype in input_control.txt.

E v

• E: Young's modulus.
• v: Poisson's ratio.

E1 E2 E3 G23 G31 G12 v23 v13 v12 angle

• E1, E2, E3: Young's moduli.
• G23, G31, G12: shear moduli.
• v23, v13, v12: Poisson's ratios.
• angle: orientation in x1-x2 plane in degrees, directed counterclockwisely from positive x1-axis to the principal material direction.

C11 C21 C31 C41 C51 C61
C12 C22 C32 C42 C52 C62
C13 C23 C33 C43 C53 C63
C14 C24 C34 C44 C54 C64
C15 C25 C35 C45 C55 C65
C16 C26 C36 C46 C56 C66

• Cij: elements of stiffness matrix.

S11 S21 S31 S41 S51 S61
S12 S22 S32 S42 S52 S62
S13 S23 S33 S43 S53 S63
S14 S24 S34 S44 S54 S64
S15 S25 S35 S45 S55 S65
S16 S26 S36 S46 S56 S66

• Sij: elements of compliance matrix.

E1 E2 G12 v12 angle

• E1, E2,: Young's moduli.
• G12: shear modulus.
• v12: Poisson's ratio.
• angle: fiber orientation in degrees, directed counterclockwisely from x1-axis to the principal material direction.

The first line defines the number of layers and materials, followed by elastic properties of each material, and the rest are layer arrangements. The first layer has the smallest $z$. The reference plane is the mid-plane.

nLayer nMat
E1 E2 G12 v12
[E1 E2 G12 v12]
[...]
mat angle thk
[mat angle thk]
[...]

E.g.

4 2
138e9 9e9 6.9e9 0.3
9e9 9e9 6e9 0.25
1 45 0.001
2 0 0.001
2 45 0.001
1 -45 0.001

In this example, there are 4 layers and 2 materials. Material 1 is defined by 138e9 9e9 6.9e9 0.3, and material 2 by 9e9 9e9 6e9 0.25. Layer 1 is defined by 1 45 0.001; it is made of material 1 with fiber angle 45° and thickness 0.001, and the range of its $z$ coordinate is [-0.002, -0.001].

• nLayer: total number of layers.
• nMat: total number of materials.
• E1, E2: Young's moduli.
• G12: shear modulus.
• v12: Poisson's ratio.
• mat: material number of this layer. E.g. mat=1 is the 1st material defined in this file.
• angle: fiber orientation in degrees, directed counterclockwisely from x1-axis to the principal material direction.
• thk: layer thickness.

Thermal conductivities and thermal expansion coefficients. Its format depends on Ttype in input_control.txt. This is only a part of valid formats. Consult §3.1.4 of Hwu (2021) for the full list.

k a

• k: thermal conductivity.

• a: thermal expansion coefficient $\alpha$.

k1 k2 k3 a1 a2 a3

• k1, k2, k3: thermal conductivities in x1, x2 and x3 directions.

• a1, a2, a3: thermal expansion coefficients $\alpha$ in x1, x2 and x3 directions.

k11 k12 k13
k12 k22 k23
k13 k23 k33
b11 b12 b13
b12 b22 b23
b13 b23 b33
[b14 b24 b34]
[b15 b25 b35]

• kij: thermal conductivities.

• bij: thermal moduli $\beta_{ij}$.

k11 k12 k13
k12 k22 k23
k13 k23 k33
a11 a12 a13
a12 a22 a23
a13 a23 a33
[a14 a24 a34]
[a15 a25 a35]

• kij: thermal conductivities.

• aij: thermal expansion coefficients $\alpha_{ij}$.

input_loadstr.txt provides the parameters used in the solutions, e.g., value of stress at infinity, location and value of point force, and parameters of elliptical hole/inclusion; therefore, its format varies for each Ltype. loadstr stands for "load and structure".

Only those I use frequently are listed here, for a full list, see Appendix E.4 of Hwu (2021).

Its format depends on Ltype in input_control.txt.

For BEM, Ltype consists of two parts: prefix (P) and fundamental solution (F).

Ltype = P + F.

For example, if P=3 and F=6, then Ltype=36. The prefix is optional; if there is no prefix and F=6, then Ltype=6.

elemType GausPts x0 y0 angle

• elemType: defines the element type. Valid labels are:

  • 1: linear element.

  • 2: quadratic element.

  • 3: linear element with cubic deflection.

• GausPts: number of Gaussain points used for the line integral.

• x0, y0: no use, but required.

• angle: no use, but required.

elemType GausPts x0 y0 angle a b

• elemType: defines the element type. Valid labels are:

  • 1: linear element.

  • 2: quadratic element.

  • 3: linear element with cubic deflection.

• GausPts: number of Gaussain points used for the line integral.

• x0, y0: center of ellipse.

• angle: angle from positive x1-axis to the major axis of the ellipse, in degrees.

• a, b: half major and minor axes of the ellipse.

elemType GausPts x0 y0 angle a b ns

• elemType: defines the element type. Valid labels are:

  • 1: linear element.

  • 2: quadratic element.

  • 3: linear element with cubic deflection.

• GausPts: number of Gaussain points used for the line integral.

• x0, y0: center of ellipse.

• angle: angle from positive x1-axis to the major axis of the ellipse, in degrees.

• a, b: half major and minor axes of the ellipse.

• ns: numer of terms for estimating infinite series, suggested value = 20.

nSub Ltype1 Ltype2 [Ltype3 ...]
loadstr1 mat1 [mat2]
loadstr2 mat1 [mat2]
[loadstr3 mat1 [mat2]]
[...]

• nSub: number of subregions

• Ltype1, Ltype2, Ltype3, ... : Ltype of each subregion.

• loadstr1, loadstr2, loadstr3, ... : input for each subregion, format depends on the corresponding Ltype. E.g. if Ltype1=4, then loadstr1 should be

  elemType GausPts x0 y0 angle a b
  

• mat1, mat2: the material for this subregion, mat2 is only required for bimaterial and inclusion. mat1=1 refers to the first material defined in input_control.txt.

loadLabel val1 [val2 val3 ...]

Avaliable loadLabel and their corresponding val are:

• Unidirectional tension with an angle, loadLabel=1

  1 sigma angle

  • sigma: value of unidirectional tension.
  • angle: angle of the load, in degrees, counterclockwise from positive x1-axis.

• Biaxial tension, loadLabel=2

  2 sigma1 sigma2

  • sigma1, sigma2: value of tension loads.

• Inplane shear, loadLabel=3

  3 sigma

  • sigma: value of inplane shear load.

• Antiplane shear, loadLabel=4

  4 sigma13 sigma23

  • sigma13, sigma23: value of anti-plane shear load $\tau_{13}$ and $\tau_{23}$.

• Stress components, loadLabel=5

  5 sigma11 sigma22 sigma23 sigma13 sigma12

  • sigma11, sigma22, sigma23, sigma13, sigma12: stress at infinity $\sigma_{11}, \sigma_{22}, \sigma_{23}, \sigma_{13}, \sigma_{12}$.

a b loadLabel val1 [val2 val3 ...]

• a, b: semi-major and minor axes of the ellipse

• loadLabel: valid labels are: 1, 2, 3, 4, 5. See Ltype 411 for loadLabels and their corresponding val1 val2 ....

a b p1 p2 p3 [p4] x1 x2

• a, b: semi-major and minor axes of the ellipse

• p1, p2, p3, p4: components of point force, p4 for piezoeletric material.

• x1, x2: (x1, x2) is the location of the point force.

a c epsilon k loadLabel val1 [val2 val3 ...]

• a, c, epsilon, k: parameters defining hole contour by $x_1 = a(\cos\psi+\epsilon\cos k\psi),\ x_2=a(c\sin\psi-\epsilon\sin k\psi).$
• loadLabel: valid labels are: 1, 2, 3, 4, 5. See Ltype 411 for loadLabels and their corresponding val1 val2 ....

Despite its mysterious name, input_variable.txt defines the points at which the results (e.g. displacement, stress, strain) are calculated and output to result.txt and figures.

Its format depends on Otype in input_control.txt. This is only a portion of all the available options. Consult Hwu (2021) for the full list.

For Otype=1,2, the first number in input_variable.txt controls the type of curve/surface. For example, if Otype=1, this number is 1 for piecewise line segments, and 2 for arc. Don't forget to include it in this file.

For multi-region problems (BFEM), Otype must be the same for each subregion. Output settings for all subregions are concatenated vertically in input_variable.txt.

Example (Otype=1)

Output setting for subregion #1:

1 2 0 0 1 1 10 3 2 10

Output setting for subregion #2:

2 0 0 5 0 360 180

"input_variable.txt":

1 2 0 0 1 1 10 3 2 10
2 0 0 5 0 360 180

• Piecewise line segments

  1 nLines x1 y1 x2 y2 nPts [x3 y3 nPts ...]

  • nLines: total number of line segments.
  • x1, y1, x2, y2: the line segment starts from (x1, y1) and ends at (x2, y2).
  • nPts: total number of points on this line segment, including its end points.

• Arc

  2 x0 y0 r start_angle end_angle nPts

  • x0, y0: center of the arc.
  • r: radius.
  • start_angle, end_angle: starting and ending angles of the arc, in degrees.
  • nPts: total number of points on the arc, end points included.

• Slanted elliptical curve

  3 x0 y0 a b start_angle end_angle nPts slant_angle

  • x0, y0: center of the ellipse
  • a, b: lengths of half major and minor axes.
  • start_angle, end_angle: range of the ellipse parameter $\psi$, in degrees.
  • slant_angle: rotation angle of the curve, counterclockwise from positive x1-axis.

• Slanted polygon-like curve

  4 x0 y0 a c epsilon k psi_0 psi_1 nPts slant_angle

  • x0, y0: center of the ellipse

  • a, c, epsilon, k: parameters of the curve

    $x_1 = a(\cos\psi+\epsilon\cos k\psi),\ x_2=a(c\sin\psi-\epsilon\sin k\psi).$

  • psi_0, psi_1: range of the ellipse parameter $\psi$, in degrees.

  • slant_angle: rotation angle of the curve, counterclockwise from positive x1-axis.

• Slanted rectangle

  1 x1 y1 x2 y2 nXPts nYPts angle

  • x1, y1, x2, y2: two corners, (x1, y1) and (x2, y2), defining the rectangle.

  • nXPts, nYPts: number of points divided in x and y directions; nXPts*nYPts points in total.

  • angle: slant angle, in degrees, counterclockwise. I guess it rotates the rectangle, defined by x1, y1, x2, and y2, with respect to its center.

• Circular Disk, ring, or sector.

  2 x0 y0 r1 theta1 r2 theta2 nRPts nThPts

  • x0, y0: center of the disk/ring.
  • r1, r2: inner and outer radii.
  • theta1, theta2: starting and ending angle, in degrees.
  • nRPts, nThPts: number of points divided in r and theta directions.

• Slanted polygon-like area

  $x_1 = a(\cos\psi+\epsilon\cos k\psi),\ x_2=a(c\sin\psi-\epsilon\sin k\psi),\ a_1 \geq a \geq a_2,\ \psi_1 \geq \psi \geq \psi_2.$

  4 x0 y0 c ep k a1 psi1 a2 psi2 naPts npsiPts angle

  • x0, y0: center of the area.
  • c, ep, k: shape parameters.
  • a1, a2: $a_1$, $a_2$, range of $a$.
  • psi1, psi2: $\psi_1$, $\psi_2$, range of $\psi$, in degrees.
  • naPts, npsiPts: number of points divided in $a$ and $\psi$ directions.
  • angle: slant angle, in degrees, counterclockwise.

nPts
x y
[x y]
[...]

• nPts: total number of points.

• x, y: coordinate of the points.

Coordinates of nodes. Arranged clockwise. Each row contains the coordinates of one node. The nodes are numbered automatically from 1 to total number of nodes.

x1 x2 [x3]
x1 x2 [x3]
[...]

Though it is called input_node1.txt, this file defines the elements by nodal connectivity. Each row defines one element.

The number 1 in the file name is its region number.

n1 n2 [n3 n4 n5 n6 n7 n8]
n1 n2 [n3 n4 n5 n6 n7 n8]
[...]

• n1: index of the first node of the element.

• n2, n3, n4, n5, n6, n7, n8: index of the second (n2) to eighth node (n8) of the element, if any.

Defines the boundary conditions of nodes. Must be defined for all the nodes.

label1 label2 label3 [label4 ...] value1 value2 value3 [value4 ...]

• label1, label2, label3, ... : type of boundary condition of the 1st to 8th degree of freedom, if any. Valid labels are:

  • 0: traction-prescribed.

  • 1: displacement-prescribed.

  • 2: contact.

• value1, value2, value3, ... : prescribed values of the 1st to 8th degree of freedom, if any.

Temperature $T$ and its derivatives $T_{,1}$ and $T_{,2}$ at each node. Each line contains data of one node in the same order as input_node1.txt.

In newer versions, modification for themalelasticity with multiple regions (BFEM) is added. Because the temperature gradients are often discontinous across material interfaces, each region needs its own temperature data.

These files should be named sequencially as input_temp1.txt, input_temp2.txt , ..., etc. All files should have exact same number of lines. The values can be set as zero for nodes outside the region.

T T1 T2
[T T1 T2]
[...]

• T: temperature.

• T1, T2: temperature gradient $T_{,1}$ and $T_{,2}$.