fem_tunnel

Outline of a program 'f90_fem_truss.f90'
Format of input data file ('csv' format)
Format of output file ('prn' format)
Sample drawings
Programs and sample data

Outline of a program 'f90_fem_ini4.f90'

This is a program for 2D-tunnel excavation analysis. Elastic and small displacement problems can be treated.
4 node isoparametric element with 4 Gauss points is used as a 4 node element. 1 node has 2 degrees of freedom in horizontal direction and vertical direction.
This program can do plane stress analysis or plane strain analysis.
In coordinate system, x-direction is defined as right-direction, y-direction is defined as upward direction and counterclockwise direction is defined as positive in rotation.
Simultaneous linear equations are solved using Cholesky method for banded matrix.
Input file name can be defined arbitrarily with the format of 'csv' and Output file format is blank selected value.
Used language for program is 'Fortran 90' and used compiler is 'GNU gfortran.'

Workable condition
Item	Description
Element	Quadrilateral element
Load	Equivalent excavation loads will be calculated automatically
Displacement	Specify the nodes and displacements at the nodes. Any values can be applied including zero
Initial stress	Initial stress in the ground shall be specified for calculation

Steps for calculation

(1) Input model caracteristics

(2) Set initial stresses

(3) Calculate equivalent excavation force

(4) Calculate displacements and stresses by loading equivalent excavation forcr

Commands for calculation

gfortran -o f90_fem_ini4 f90_fem_ini4.f90
./f90_fem_ini4  gamma  sig0  r  fnameR.csv  fnameW.prn

f90_fem_ini4	Program for execution
gamma	Unit weight of ground material
sig0	Initial vertical stress at y=0 on the model
r	Lateral pressure ratio (lateral stress / vertical stress)
fnameR.csv	Input file (csv format)
fnameW.prn	Output file (prn format: space separated values)

Format of input data file ('csv' format)

Comment                               # Comments
NODT,NELT,MATEL,KOX,KOY,NF,NSTRES,IPR # Basic values for analysis
t,Em,po,c,phi,sigt                    # material properties
    ....(1 to MATEL)....              #
node-1,node-2,node-3,node-4,matno     # Element connectivity, material set number
    ....(1 to NELT)....               #
x,y                                   # Node coordinates
    ....(1 to NODT)....               #
nokx,rdisx                            # Restricted node number and displacement in x-direction
    ....(1 to KOX)....                # (Omit data input if KOX=0)
noky,rdisy                            # Restricted node number and displacement in y-direction
    ....(1 to KOY)....                # (Omit data input if KOY=0)
nonf                                  # loaded node number
    ....(1 to NF)....                 # (Omit data input if NF=0)

nod	: Number of nodes of each element (3 or 4)	t	: Thickness of element
NODT	: Number of nodes	Em	: Elastic modulus of element
NELT	: Number of elements	po	: Poisson's ratio of element
MATEL	: Number of material set	c	: Pure shearing strength of element
KOX	: Number of restricted nodes in x-direction	phi	: Internal friction angle of element
KOY	: Number of restricted nodes in y-direction	sigt	: Tensile strength of element
NF	: Number of loaded nodes	matno	: Material set number
NSTRES	: 0: plane strain, 1: plane stress
IPR	: Output format of stresses (0 or 1)
	(0: All Gauss points, 1: average)

Notice

Restricted node means the node which has known (given) displacement. As a known (given) value of nodal displacement, any value can be given including zero for a restricted node.

Format of output file ('prn' format)

Comment
nod  NODT  NELT  MATEL  KOX  KOY  NF  NSTRES  IPR
(Each value for above item)
*node characteristics
node  x  y  fx  fy  fix-x  fix-y  rdis-x  rdis-y
    node   : Number of nodes
    x,y    : Coordinates of x & y-direction
    fx,fy  : Loads of x & y-direction
    fix-x  : x-direction restricted condition (1: restricted, 0: not restricted)
    fix-y  : y-direction restricted condition (1: restricted, 0: not restricted)
    rdis-x : x-direction specified displacement
    rdis-y : y-direction specified displacement
  .....(1 to NODT).....
*element characteristics
element  node-1  node-2  node-3  node-4  E  po  t  chs  phi  sigt  matno
    element : Element number
    node-1,node-2,node-3,node-4 : Element-nodes relationship
    E       : Elastic modulus of element
    po      : Poisson's ratio of element
    t       : Thickness of element
    chs     : Pure shearing strength of element
    phi     : Internal friction angle of element
    sigt    : Tensile stress of element
    matno   : material set number
  .....(1 to NELT).....
*displacements and forces
node  coord-x  coord-y  dis-x  dis-y  reac-x  reac-y  ftvec-x  ftvec-y
    node            : Node number
    coord-x,coord,y : x & y-coordinates
    dist-x,dist-y   : x & y-direction displacements
    reac-x,reac-y   : These values shall be zero.
    ftvec-x,ftvec-y : x & y-direction loaded forces
  .....(1 to NODT).....
*stresses
element  kk  coord-x  coord-y  sig-x  sig-y  tau-xy  ps1  ps2  ang  SF  matno
    element            : Element number
    kk                 : Gauss point number (IPR=0: 1 to 4, IPR=1: Average stress of element)
    coord-x,coord-y    : Coordinates of calculation point of stresses
    sig-x,sig-y,tau-xy : Normal stresses in x & y-direction, shearing stress
    ps1,ps2,ang        : 1st & 2nd principal stresses, principal direction (degree)
    SF                 : Point safety factor
    matno              : Material set number
  .....(1 to NELT, depending on the value of IPR).....
#,Summary
#,NELT=(Number of elements) NODT=(Number of nodes) nt=(nt) mm=(mm) ib=(ib)
#,ps1_range, (ps1_min), (ps1_max)
#,ps2_range, (ps2_min), (ps2_max)
#,sf_range,  (sf_min) , (sf_max)
#,dis_max,  (dis_max)
    nt : Total degrees of freedom of FE equation
    mm : Dimension of reduced FE equation
    ib : band width of reduced FE equation

Sample drawings

Mesh by GMT

fig_gmt_mesh.png	fig_gmt_disp.png

Contours by GMT (unit of stresses: tf/m$^2$)

fig_gmt_ps1.png	fig_gmt_ps2.png	fig_gmt_sf.png

Arrow diagrams by GMT

fig_gmt_vec1.png	fig_gmt_vec2.png

Contours by Python3

fig_ini4_cont1.png	fig_ini4_cont2.png	fig_ini4_cont3.png

Arrow diagrams by Python3

fig_ini4_vect1.png	fig_ini4_vect2.png

Programs and sample data

Fortran programs

Filename	Description
f90_fem_ini4.f90	Tunnel Excavation Analysis
f90_gmt_ini4_con.f90	Drawing using GMT (principal stress and safety factor contour)
f90_gmt_ini4_vec.f90	Drawing using GMT (stress vector arrow diagram)
py_fem_ini4_cont.py	Drawing using Python3 (principal stress and safety factor contours)
py_fem_ini4_vect.py	Drawing using Python3 (stress vector arrow diagram)

Usage of Fortran drawing programs

./f90_gmt_ini4_con fnameW.prn
./f90_gmt_ini4_vec fnameW.prn

where, fnameW.prn is output data file by f90_fem_ini4.

F90 program	Output drawing (fixed name)	Description
f90_gmt_ini4_con.f90	_fig_gmt_mesh.eps	Mesh drawing
	_fig_gmt_ps1.eps	First principal stress contour
	_fig_gmt_ps2.eps	Second principal stress contour
	_fig_gmt_sf.eps	Safety factor contour
f90_gmt_ini4_vec.f90	_fig_gmt_vec1.eps	First principal stress arrow diagram
	_fig_gmt_vec2.eps	Second principal stress arrow diagram
	_fig_gmt_disp.eps	Displacement mode

Usage of Python3 drawing programs

python3 py_fem_ini4_cont.py fnameW.prn
python3 py_fem_ini4_vect.py fnameW.prn