# Dynamic Crack Analysis in Isotropic/Orthotropic Media via Extended Isogeometric Analysis

Authors: Tiantang Yu; Yongling Lai; Shuohui Yin
Journal: Mathematical Problems in Engineering
 (2014)

---

## Abstract

The extended isogeometric analysis (X-IGA) is the combination of the extended finite element method (X-FEM) and the isogeometric analysis (IGA), so the X-IGA possesses the advantages of both methods. In this paper, the X-IGA is extended to investigate the dynamic stress intensity factors of cracked isotropic/orthotropic media under impact loading. For this purpose, a corresponding dynamic X-IGA model is developed, the Newmark time integration scheme is used to achieve a dynamic response, and the dynamic stress intensity factors are evaluated through the contour interaction integral technique. Numerical simulations show that the X-IGA results agree with other available reference solutions, and accurate results can be obtained by using the X-IGA with a relatively coarse mesh.


---

## Body

## 1. Introduction

Among existing numerical methods, the extended finite element method (X-FEM), which was developed in 1999 by Belytschko and his coworkers [1, 2], is the most effective method for solving discontinuity problems. In the X-FEM, the geometry of the discontinuity is independent of the computational mesh; therefore, the computational mesh does not need to be updated in the simulation of the propagation of discontinuities. In the past decades, the extensive literature on improving and applying the original X-FEM in modeling discontinuities has been published [3–15]. The X-FEM aims to enrich standard finite element approximation by using discontinuous basis functions in the framework of partition of unity.The standard finite element approximation is element-based polynomial approximation, and it yields discretization errors in complex geometry. In order to overcome the drawback, the isogeometric analysis (IGA) was proposed by Hughes et al. [16] in 2005. The principle of the IGA is that nonuniform rational B-splines (NURBS) basis functions are employed as shape functions for geometric description and field approximation. The IGA has some unique advantages [16], so it has been developed and applied in many fields including structural mechanics [17–24], solid mechanics [25–27], fluid mechanics [28], and contact mechanics [29].Recently, the IGA has been improved with enrichment functions to solve linear elastic fracture mechanics problems [30–32] and curved interface problems [33]. The extended isogeometric analysis (X-IGA) was proposed in [31] and contains the inherent advantages of both IGA and X-FEM. Currently, the X-IGA is further developed to analyze cracked orthotropic media [34]. In this study, the X-IGA is extended to investigate the dynamic fracture behavior of stationary cracks in isotropic/orthotropic media. A corresponding dynamic X-IGA model is developed, and the Newmark time integration scheme is used to achieve a dynamic response. The dynamic stress intensity factors (DSIFs) are evaluated by using the contour interaction integral technique.The paper is organized as follows. Section2 briefly reviews the fundamentals of NURBS-based IGA. The X-IGA for dynamic cracked isotropic/orthotropic media is described in Section 3. The DSIFs are derived by using the contour interaction integral technique in Section 4. Section 5 presents the numerical results obtained via X-IGA and compares the results with other solutions. Finally, conclusions and prospects are drawn in Section 6.


## 2. NURBS-Based Isogeometric Analysis

In the CAD, a two-dimensional NURBS surface can be constructed as [35]
(1)
S
(
ξ
,
η
)
=
∑
i
=
1
n
∑
j
=
1
m
R
i
,
j
p
,
q
(
ξ
,
η
)
B
i
,
j
,
where R
i
,
j
p
,
q
(
ξ
,
η
) is the NURBS basis function; p and q are, respectively, the orders of basis functions in the ξ direction and the η direction; B
i
,
j represents the coordinate of control point; and n and m are the numbers of basis functions in the ξ direction and the η direction, respectively.In NURBS-based isogeometric analysis [16], the displacement field u is approximated similarly to the isoparametric finite element method as follows:
(2)
u
=
∑
I
=
1
n
R
I
(
ξ
)
u
I
,
where n is the number of the control points; ξ is the parametric coordinate; R
I
(
ξ
) denotes the NURBS basis function at control point I; and u
I is the displacements of control point I.


## 3. X-IGA for Dynamic Cracked Media

### 3.1. Displacement Approximation

Similar to the X-FEM, the X-IGA aims to enrich the standard IGA approximation by using additional functions on the basis of the partition of unity to model discontinuities. For crack problems, the enriched displacement approximation can be expressed as(3)
u
h
(
x
)
=
∑
i
∈
N
s
R
i
(
x
)
u
i
+
∑
j
∈
N
cut
R
j
(
x
)
[
H
(
x
)
-
H
(
x
j
)
]
a
j
+
∑
k
∈
N
tip
R
k
(
x
)
∑
α
=
1
4
[
F
α
(
x
)
-
F
α
(
x
k
)
]
b
k
α
,
where R
i
(
x
), R
j
(
x
), and R
k
(
x
) are the shape functions of standard IGA (NURBS basis functions); u
i, a
j, and b
k
α are, respectively, the displacement and enrichment variable vectors at control point; N
s, N
cut, and N
tip are, respectively, the set of all control points in the computational domain, the set of control points enriched with a modified Heaviside step function H
(
x
), and the set of control points enriched with the crack-tip branch enrichment functions F
α
(
x
); (
α
=
1
,
…
,
4
); and the basis function support of the control point in N
cut is completely split by the crack, whereas that of the control point in N
tip is partly split by the crack.The modified Heaviside step functionH
(
x
) is given by
(4)
H
(
x
)
=
{
+
1
x
>
0
,
-
1
x
<
0
.The crack-tip branch enrichment functionsF
α
(
x
) for orthotropic materials are defined as [36](5)
{
F
α
(
x
)
}
α
=
1
4
=
{
r
cos
⁡
θ
1
2
g
1
(
θ
)
r
cos
⁡
θ
2
2
g
2
(
θ
)
r
sin
θ
1
2
g
1
(
θ
)
r
sin
θ
2
2
g
2
(
θ
)
}
,where r and θ are the crack-tip local polar coordinates, and
(6)
g
k
(
θ
)
=
(
cos
⁡
θ
+
μ
k
x
sin
θ
)
2
+
(
μ
k
y
sin
θ
)
2
,
hhl
θ
k
(
θ
)
=
arctg
(
μ
k
y
sin
θ
cos
⁡
θ
+
μ
k
x
sin
θ
)
,
where μ
k
x and μ
k
y
(
k
=
1,2
) are crack-tip material parameters and μ
k
=
μ
k
x
+
i
μ
k
y  (
k
=
1,2
) are the roots of
(7)
a
11
μ
4
-
2
a
16
μ
3
+
(
2
a
12
+
a
66
)
μ
2
-
2
a
26
μ
+
a
22
=
0
,
where a
i
j is the compliance coefficients.Equation (5) cannot be directly used for isotropic materials because of the presence of 0/0 in the equation. The crack-tip branch enrichment functions F
α
(
x
) for isotropic materials are defined as [1]
(8)
{
F
α
(
x
)
}
α
=
1
4
=
{
r
sin
θ
2
r
cos
⁡
θ
2
r
sin
θ
2
sin
θ
r
cos
⁡
θ
2
sin
θ
}
.


### 3.2. Discrete Equilibrium Equations

Considering a bodyΩ with an initial traction-free crack in the state of dynamic equilibrium, the weak form of the momentum equation is
(9)
∫
Ω
ρ
u
¨
·
δ
u
d
Ω
+
∫
Ω
σ
·
δ
ε
d
Ω
=
∫
Ω
f
·
δ
u
d
Ω
+
∫
Γ
f
¯
·
δ
u
d
Γ
,
where u is the displacement vector; σ and ε are the stress and strain tensors, respectively; f and f
¯ are the body force and external traction vectors, respectively; and ρ is the mass density.The discretized form of (9) using the X-IGA approximation (3) can be written as
(10)
M
δ
¨
+
K
δ
=
F
,
where δ
=
[
u
a
b
]
T is the vector of the unknown variable at control points and M, K, and F are, respectively, the global mass matrix, stiffness matrix, and external force vector at control points.The element contribution toM is
(11)
m
i
j
=
[
m
i
j
u
u
m
i
j
u
a
m
i
j
u
b
m
i
j
a
u
m
i
j
a
a
m
i
j
a
b
m
i
j
b
u
m
i
j
b
a
m
i
j
b
b
]
,
where
(12)
m
i
j
u
u
=
∫
Ω
ρ
(
R
i
)
T
R
j
d
Ω
,
m
i
j
a
a
=
∫
Ω
ρ
(
H
(
x
)
-
H
(
x
i
)
)
(
H
(
x
)
-
H
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
,
m
i
j
b
b
=
∫
Ω
ρ
(
F
α
(
x
)
-
F
α
(
x
i
)
)
(
F
α
(
x
)
-
F
α
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
,
hhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh
(
α
=
1
~
4
)
,
m
i
j
u
a
=
∫
Ω
ρ
(
H
(
x
)
-
H
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
,
m
i
j
u
b
=
∫
Ω
ρ
(
F
α
(
x
)
-
F
α
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
(
α
=
1
~
4
)
,
m
i
j
a
b
=
∫
Ω
ρ
(
H
(
x
)
-
H
(
x
i
)
)
(
F
α
(
x
)
-
F
α
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
hhhhlhhhhhhhhhhhhhhhhhhhhhhhhhhhhh
(
α
=
1
~
4
)
,
and the element contribution to K is
(13)
k
i
j
=
[
k
i
j
u
u
k
i
j
u
a
k
i
j
u
b
k
i
j
a
u
k
i
j
a
a
k
i
j
a
b
k
i
j
b
u
k
i
j
b
a
k
i
j
b
b
]
,
where
(14)
k
i
j
r
s
=
∫
Ω
e
(
B
i
r
)
T
D
B
j
s
d
Ω
(
r
,
s
=
u
,
a
,
b
)
with
(15)
B
i
u
=
[
R
i
,
x
0
0
R
i
,
y
R
i
,
y
R
i
,
x
]
,
B
i
a
=
(
H
(
x
)
-
H
(
x
i
)
)
[
R
i
,
x
0
0
R
i
,
y
R
i
,
y
R
i
,
x
]
,
B
i
b
=
[
B
i
b
1
B
i
b
2
B
i
b
3
B
i
b
4
]
,
B
i
b
α
=
[
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
x
0
0
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
y
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
y
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
x
]
hhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh
(
α
=
1
~
4
)
.Additionally, the element contribution toF is
(16)
f
i
=
[
f
i
u
f
i
a
f
i
b
]
with
(17)
f
i
u
=
∫
∂
Ω
e
R
i
t
¯
d
Γ
+
∫
Ω
e
R
i
t
d
Ω
,
f
i
a
=
∫
∂
Ω
e
R
i
(
H
(
x
)
-
H
(
x
i
)
)
t
¯
d
Γ
+
∫
Ω
e
R
i
(
H
(
x
)
-
H
(
x
i
)
)
t
d
Ω
,
f
i
b
α
=
∫
∂
Ω
e
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
t
¯
d
Γ
+
∫
Ω
e
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
t
d
Ω
(
α
=
1
~
4
)
,
where t
¯ and t are the external traction and body force vectors, respectively.


### 3.3. Numerical Integration Scheme

The Gauss quadrature scheme is employed in the X-IGA. To obtain an accurate integration for crack tip elements and elements cut by crack, the triangular subdomain technique is used in the same way as that of the X-FEM [1]. For additional details, refer to [31].


### 3.4. Time Integration Scheme

The Newmark method is adopted for the time integration of (10). At time step n, the discrete simultaneous equations are described as
(18)
(
M
+
β
Δ
t
2
K
)
δ
¨
n
=
F
n
-
K
[
Δ
t
2
2
δ
n
-
1
+
Δ
t
δ
˙
n
-
1
hhhhhhhhhhhhhhhhhhhhh
+
(
1
-
2
β
)
Δ
t
2
2
δ
¨
n
-
1
]
,
hhhhhhhhhh
δ
˙
n
=
δ
˙
n
-
1
+
[
(
1
-
α
)
δ
¨
n
-
1
+
α
δ
¨
n
]
Δ
t
,
hhhhhhhhhh
δ
n
=
δ
n
-
1
+
δ
˙
n
-
1
Δ
t
hhhhhhhhhhhhhh
+
[
(
1
2
-
β
)
δ
¨
n
-
1
+
β
δ
¨
n
]
Δ
t
2
,
where Δ
t is the time step and the unconditionally stable conditions are α
≥
0.5 and β
≥
0.25
(
α
+
0.5
)
2. In this study, α
=
0.5 and β
=
0.25 are chosen.The implementation procedure of the Newmark time integration scheme applying into the X-IGA is outlined as follows.(1)
The initial velocities and initial displacements are set to be zero; then the initial accelerations are obtained with (10).
(2)
Compute the constantsa
0
=
1
/
β
Δ
t
2, a
1
=
α
/
β
Δ
t, a
2
=
1
/
β
Δ
t, a
3
=
(
1
/
2
β
)
-
1, a
4
=
(
α
/
β
)
-
1, a
5
=
(
Δ
t
/
2
)
(
(
α
/
β
)
-
2
), a
6
=
Δ
t
(
1
-
α
), and a
7
=
α
Δ
t.
(3)
ComputeK
¯
=
a
0
M
+
K.
(4)
Loop over the time steps as follows:
(a)
computeF
n
¯
=
F
n
+
M
(
a
0
δ
n
-
1
+
a
2
δ
˙
n
-
1
+
a
3
δ
¨
n
-
1
);
(b)
compute the displacementsδ
n by solving K
¯
δ
n
=
F
¯
n;
(c)
solveδ
¨
n
=
a
0
(
δ
n
-
δ
n
-
1
)
-
a
2
δ
˙
n
-
1
-
a
3
δ
¨
n
-
1 for the value δ
¨
n, and the velocities δ
˙
n are calculated by using δ
˙
n
=
δ
˙
n
-
1
+
a
6
δ
¨
n
-
1
+
a
7
δ
¨
n.
(5)
Repeat the loop for the next time step until the maximum time step is reached.


## 3.1. Displacement Approximation

Similar to the X-FEM, the X-IGA aims to enrich the standard IGA approximation by using additional functions on the basis of the partition of unity to model discontinuities. For crack problems, the enriched displacement approximation can be expressed as(3)
u
h
(
x
)
=
∑
i
∈
N
s
R
i
(
x
)
u
i
+
∑
j
∈
N
cut
R
j
(
x
)
[
H
(
x
)
-
H
(
x
j
)
]
a
j
+
∑
k
∈
N
tip
R
k
(
x
)
∑
α
=
1
4
[
F
α
(
x
)
-
F
α
(
x
k
)
]
b
k
α
,
where R
i
(
x
), R
j
(
x
), and R
k
(
x
) are the shape functions of standard IGA (NURBS basis functions); u
i, a
j, and b
k
α are, respectively, the displacement and enrichment variable vectors at control point; N
s, N
cut, and N
tip are, respectively, the set of all control points in the computational domain, the set of control points enriched with a modified Heaviside step function H
(
x
), and the set of control points enriched with the crack-tip branch enrichment functions F
α
(
x
); (
α
=
1
,
…
,
4
); and the basis function support of the control point in N
cut is completely split by the crack, whereas that of the control point in N
tip is partly split by the crack.The modified Heaviside step functionH
(
x
) is given by
(4)
H
(
x
)
=
{
+
1
x
>
0
,
-
1
x
<
0
.The crack-tip branch enrichment functionsF
α
(
x
) for orthotropic materials are defined as [36](5)
{
F
α
(
x
)
}
α
=
1
4
=
{
r
cos
⁡
θ
1
2
g
1
(
θ
)
r
cos
⁡
θ
2
2
g
2
(
θ
)
r
sin
θ
1
2
g
1
(
θ
)
r
sin
θ
2
2
g
2
(
θ
)
}
,where r and θ are the crack-tip local polar coordinates, and
(6)
g
k
(
θ
)
=
(
cos
⁡
θ
+
μ
k
x
sin
θ
)
2
+
(
μ
k
y
sin
θ
)
2
,
hhl
θ
k
(
θ
)
=
arctg
(
μ
k
y
sin
θ
cos
⁡
θ
+
μ
k
x
sin
θ
)
,
where μ
k
x and μ
k
y
(
k
=
1,2
) are crack-tip material parameters and μ
k
=
μ
k
x
+
i
μ
k
y  (
k
=
1,2
) are the roots of
(7)
a
11
μ
4
-
2
a
16
μ
3
+
(
2
a
12
+
a
66
)
μ
2
-
2
a
26
μ
+
a
22
=
0
,
where a
i
j is the compliance coefficients.Equation (5) cannot be directly used for isotropic materials because of the presence of 0/0 in the equation. The crack-tip branch enrichment functions F
α
(
x
) for isotropic materials are defined as [1]
(8)
{
F
α
(
x
)
}
α
=
1
4
=
{
r
sin
θ
2
r
cos
⁡
θ
2
r
sin
θ
2
sin
θ
r
cos
⁡
θ
2
sin
θ
}
.


## 3.2. Discrete Equilibrium Equations

Considering a bodyΩ with an initial traction-free crack in the state of dynamic equilibrium, the weak form of the momentum equation is
(9)
∫
Ω
ρ
u
¨
·
δ
u
d
Ω
+
∫
Ω
σ
·
δ
ε
d
Ω
=
∫
Ω
f
·
δ
u
d
Ω
+
∫
Γ
f
¯
·
δ
u
d
Γ
,
where u is the displacement vector; σ and ε are the stress and strain tensors, respectively; f and f
¯ are the body force and external traction vectors, respectively; and ρ is the mass density.The discretized form of (9) using the X-IGA approximation (3) can be written as
(10)
M
δ
¨
+
K
δ
=
F
,
where δ
=
[
u
a
b
]
T is the vector of the unknown variable at control points and M, K, and F are, respectively, the global mass matrix, stiffness matrix, and external force vector at control points.The element contribution toM is
(11)
m
i
j
=
[
m
i
j
u
u
m
i
j
u
a
m
i
j
u
b
m
i
j
a
u
m
i
j
a
a
m
i
j
a
b
m
i
j
b
u
m
i
j
b
a
m
i
j
b
b
]
,
where
(12)
m
i
j
u
u
=
∫
Ω
ρ
(
R
i
)
T
R
j
d
Ω
,
m
i
j
a
a
=
∫
Ω
ρ
(
H
(
x
)
-
H
(
x
i
)
)
(
H
(
x
)
-
H
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
,
m
i
j
b
b
=
∫
Ω
ρ
(
F
α
(
x
)
-
F
α
(
x
i
)
)
(
F
α
(
x
)
-
F
α
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
,
hhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh
(
α
=
1
~
4
)
,
m
i
j
u
a
=
∫
Ω
ρ
(
H
(
x
)
-
H
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
,
m
i
j
u
b
=
∫
Ω
ρ
(
F
α
(
x
)
-
F
α
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
(
α
=
1
~
4
)
,
m
i
j
a
b
=
∫
Ω
ρ
(
H
(
x
)
-
H
(
x
i
)
)
(
F
α
(
x
)
-
F
α
(
x
j
)
)
(
R
i
)
T
R
j
d
Ω
hhhhlhhhhhhhhhhhhhhhhhhhhhhhhhhhhh
(
α
=
1
~
4
)
,
and the element contribution to K is
(13)
k
i
j
=
[
k
i
j
u
u
k
i
j
u
a
k
i
j
u
b
k
i
j
a
u
k
i
j
a
a
k
i
j
a
b
k
i
j
b
u
k
i
j
b
a
k
i
j
b
b
]
,
where
(14)
k
i
j
r
s
=
∫
Ω
e
(
B
i
r
)
T
D
B
j
s
d
Ω
(
r
,
s
=
u
,
a
,
b
)
with
(15)
B
i
u
=
[
R
i
,
x
0
0
R
i
,
y
R
i
,
y
R
i
,
x
]
,
B
i
a
=
(
H
(
x
)
-
H
(
x
i
)
)
[
R
i
,
x
0
0
R
i
,
y
R
i
,
y
R
i
,
x
]
,
B
i
b
=
[
B
i
b
1
B
i
b
2
B
i
b
3
B
i
b
4
]
,
B
i
b
α
=
[
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
x
0
0
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
y
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
y
(
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
)
,
x
]
hhhhhhhhhhhhhhhhhhhhhhhhhhhhhhhh
(
α
=
1
~
4
)
.Additionally, the element contribution toF is
(16)
f
i
=
[
f
i
u
f
i
a
f
i
b
]
with
(17)
f
i
u
=
∫
∂
Ω
e
R
i
t
¯
d
Γ
+
∫
Ω
e
R
i
t
d
Ω
,
f
i
a
=
∫
∂
Ω
e
R
i
(
H
(
x
)
-
H
(
x
i
)
)
t
¯
d
Γ
+
∫
Ω
e
R
i
(
H
(
x
)
-
H
(
x
i
)
)
t
d
Ω
,
f
i
b
α
=
∫
∂
Ω
e
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
t
¯
d
Γ
+
∫
Ω
e
R
i
(
F
α
(
x
)
-
F
α
(
x
i
)
)
t
d
Ω
(
α
=
1
~
4
)
,
where t
¯ and t are the external traction and body force vectors, respectively.


## 3.3. Numerical Integration Scheme

The Gauss quadrature scheme is employed in the X-IGA. To obtain an accurate integration for crack tip elements and elements cut by crack, the triangular subdomain technique is used in the same way as that of the X-FEM [1]. For additional details, refer to [31].


## 3.4. Time Integration Scheme

The Newmark method is adopted for the time integration of (10). At time step n, the discrete simultaneous equations are described as
(18)
(
M
+
β
Δ
t
2
K
)
δ
¨
n
=
F
n
-
K
[
Δ
t
2
2
δ
n
-
1
+
Δ
t
δ
˙
n
-
1
hhhhhhhhhhhhhhhhhhhhh
+
(
1
-
2
β
)
Δ
t
2
2
δ
¨
n
-
1
]
,
hhhhhhhhhh
δ
˙
n
=
δ
˙
n
-
1
+
[
(
1
-
α
)
δ
¨
n
-
1
+
α
δ
¨
n
]
Δ
t
,
hhhhhhhhhh
δ
n
=
δ
n
-
1
+
δ
˙
n
-
1
Δ
t
hhhhhhhhhhhhhh
+
[
(
1
2
-
β
)
δ
¨
n
-
1
+
β
δ
¨
n
]
Δ
t
2
,
where Δ
t is the time step and the unconditionally stable conditions are α
≥
0.5 and β
≥
0.25
(
α
+
0.5
)
2. In this study, α
=
0.5 and β
=
0.25 are chosen.The implementation procedure of the Newmark time integration scheme applying into the X-IGA is outlined as follows.(1)
The initial velocities and initial displacements are set to be zero; then the initial accelerations are obtained with (10).
(2)
Compute the constantsa
0
=
1
/
β
Δ
t
2, a
1
=
α
/
β
Δ
t, a
2
=
1
/
β
Δ
t, a
3
=
(
1
/
2
β
)
-
1, a
4
=
(
α
/
β
)
-
1, a
5
=
(
Δ
t
/
2
)
(
(
α
/
β
)
-
2
), a
6
=
Δ
t
(
1
-
α
), and a
7
=
α
Δ
t.
(3)
ComputeK
¯
=
a
0
M
+
K.
(4)
Loop over the time steps as follows:
(a)
computeF
n
¯
=
F
n
+
M
(
a
0
δ
n
-
1
+
a
2
δ
˙
n
-
1
+
a
3
δ
¨
n
-
1
);
(b)
compute the displacementsδ
n by solving K
¯
δ
n
=
F
¯
n;
(c)
solveδ
¨
n
=
a
0
(
δ
n
-
δ
n
-
1
)
-
a
2
δ
˙
n
-
1
-
a
3
δ
¨
n
-
1 for the value δ
¨
n, and the velocities δ
˙
n are calculated by using δ
˙
n
=
δ
˙
n
-
1
+
a
6
δ
¨
n
-
1
+
a
7
δ
¨
n.
(5)
Repeat the loop for the next time step until the maximum time step is reached.


## 4. DSIF Calculations

The DSIFs are evaluated by using the domain form of the contour interaction integral technique. The following states are considered: state 1 (σ
i
j
(
1
), ε
i
j
(
1
), u
i
(
1
)), which corresponds to the actual state; state 2 (σ
i
j
(
2
), ε
i
j
(
2
), u
i
(
2
)), which is an auxiliary state that will be selected as the asymptotic field for model I or II. The interaction integral may be written as [37]
(19)
I
(
1,2
)
=
∫
A
[
-
W
(
1,2
)
δ
1
j
+
σ
i
j
(
1
)
∂
u
i
(
2
)
∂
x
1
+
σ
i
j
(
2
)
∂
u
i
(
1
)
∂
x
1
]
∂
q
∂
x
j
d
A
+
∫
A
ρ
u
¨
(
1
)
∂
u
i
(
2
)
∂
x
1
d
A
,
where W
(
1,2
)
=
σ
i
j
(
1
)
ε
i
j
(
2
)
=
σ
i
j
(
2
)
ε
i
j
(
1
) is the interaction strain energy and q is a weighing function. The second term of the right-hand side of the equation is the contribution of the inertia forces.For isotropic materials, the relation between the interaction integralI
(
1,2
) and the stress intensity factors is expressed as [37]
(20)
I
(
1,2
)
=
2
(
K
I
(
1
)
K
I
(
2
)
+
K
I
I
(
1
)
K
I
I
(
2
)
)
E
*
,
where E
*
=
E is for plane stress and E
*
=
E
/
(
1
-
ν
2
) is for plane strain.The relation between the interaction integralI
(
1,2
) and the stress intensity factors for orthotropic materials is expressed as [36]
(21)
I
(
1,2
)
=
2
c
11
K
I
(
1
)
K
I
(
2
)
+
c
12
(
K
I
(
1
)
K
I
I
(
2
)
+
K
I
(
2
)
K
I
I
(
1
)
)
+
2
c
22
K
I
I
(
1
)
K
I
I
(
2
)
with
(22)
c
11
=
-
a
22
2
Im
⁡
(
μ
1
+
μ
2
μ
1
μ
2
)
,
c
12
=
-
a
22
2
Im
⁡
(
1
μ
1
μ
2
)
+
a
11
2
Im
⁡
(
μ
1
μ
2
)
,
c
22
=
a
11
2
Im
⁡
(
μ
1
+
μ
2
)
.Choosing state 2 as mode I or mode II leads to mode-I or mode-II DSIF in terms of the interaction integral.For isotropic materials,(23)
{
K
I
(
1
)
K
I
I
(
1
)
}
=
E
*
2
{
I
(
1
,
I
)
I
(
1
,
I
I
)
}
.For orthotropic materials,(24)
{
K
I
(
1
)
K
I
I
(
1
)
}
=
[
2
c
11
c
12
c
12
2
c
22
]
-
1
{
I
(
1
,
I
)
I
(
1
,
I
I
)
}
.


## 5. Numerical Simulations

In this section, several examples of stationary cracks in isotropic/orthotropic media with available reference solutions are investigated to assess the accuracy of the proposed approach. The first two examples demonstrate the efficiency of the X-IGA for isotropic material problems, and the next two illustrate the application of the X-IGA to orthotropic material problems, whereas the last example demonstrates the capability and versatility of the X-IGA in modeling the complicated geometries. In all examples, the Heaviside step loading is considered and degree 3 X-IGA is adopted; that is,p
=
q
=
3. Meshes are generated with linear parameterization method. In the case of investigating the effects of different meshes on the results, different uniform meshes are adopted; in other cases, the fine meshes are used around the crack, while the coarse meshes are used in other domains; thus excellent accuracy can be achieved at a low cost.


### 5.1. Semi-Infinite Crack in an Infinite Isotropic Plate

For the first example, we consider an infinite isotropic plate with a semi-infinite crack subjected to a tensile stress wave perpendicular to the crack face, as shown in Figure1. The material properties of the media are as follows: Young’s modulus E
= 210 GPa, Poisson’s ratio ν
= 0.3, and mass density ρ
= 8000 kg/m3. The tensile stress σ
0
= 500 MPa. The theoretical solution of mode-I DSIF for the stationary crack is expressed as [38]
(25)
K
I
dyn
(
t
)
=
{
0
,
t
<
t
c
2
σ
0
1
-
v
c
d
(
t
-
t
c
)
(
1
-
2
v
)
π
,
t
c
≤
t
≤
3
t
c
,
where t
c
=
(
H
/
c
d
)  (c
d: the dilatational wave speed).Figure 1
Geometry and loading of a semi-infinite crack.


#### 5.1.1. Convergence Study of the DSIF versus Meshes

In this study, five meshes with 27 × 9, 47 × 17, 67 × 25, 87 × 33, and 107 × 41 elements are considered. Figure2 illustrates the typical regular mesh of 27 × 9 elements.Figure 2
Control points and elements for degree 3 X-IGA with 27 × 9 elements (the blue cross represents the control points enriched by Heaviside function; the blue square represents the control points enriched by crack tip enrichment functions; the red line represents the crack).Figure3 presents the normalized mode-I DSIFs as a function of time for the five considered meshes. In the mesh refinement process, the mode-I DSIFs obtained by the X-IGA converge well to the theoretical solution. The maximum error increases beside t
=
t
c for the coarse mesh and decreases significantly for the fine mesh. Similar observation was noted in the simulation of the singular edge-based smoothed finite element method (sES-FEM) [37]. The same simulations were conducted by other researches by using different methods. However, the same conclusion is recorded; that is, the error is lower in later stages than in earlier stages. The relative errors during 1.5
t
c
≤
t
≤
3
t
c are presented in Figure 4, in which the error is unchanged as the mesh is refined to a certain extent. In this study, the errors are the same for meshes with 67 × 25, 87 × 33, and 107 × 41 elements.Figure 3
Convergence of the normalized mode-I DSIFs versus mesh densities.Figure 4
Convergence in percentage errors for the normalized mode-I DSIFs versus mesh densities during1.5
t
c
≤
t
≤
3
t
c.


#### 5.1.2. Comparison with X-FEM Results

The X-IGA results using 67 × 25 mesh elements are compared with three available X-FEM solutions, including the X-FEM using mass lumping (120 × 60 quadrilateral elements) [39], the mesh-free enriched X-FEM with crack-tip enrichments (78 × 39 quadrilateral elements) [38], and X-FEM with a new enrichment function (120 × 59 quadrilateral elements) [40]. Figure 5 shows the comparative study of the normalized mode-I DSIFs. The results obtained with X-FEM using mass lumping [39] are unstable, and the accuracy of this method is the worst recorded among the four approaches. The results obtained with the X-IGA, the X-FEM using an enrichment function, and the mesh-free enriched X-FEM with crack-tip enrichments are very accurate during 1.5
t
c
≤
t
≤
3
t
c; however, the results are relatively poor during t
≤
1.5
t
c. The X-IGA results are more accurate than other numerical results during 1.5
t
c
≤
t
≤
3
t
c.Figure 5
Comparison of the normalized mode-I DSIFs among different methods.


### 5.2. Arbitrarily Oriented Central Crack in a Rectangular Isotropic Plate

We consider a rectangular isotropic plate with an arbitrarily oriented central crack subjected on the top and bottom of a uniform impact loading, as shown in Figure6(a). The length of the crack is 2
a
= 4.8 mm. The material properties of the plate are as follows: Young’s modulus E
= 200 GPa, Poisson’s ratio ν
= 0.3, and mass density ρ
= 5000 kg/m3. Three crack inclination angles are investigated, and the total time of the simulation is 20 μs. In the X-IGA, the computational mesh is independent of the crack; therefore, the computation meshes for the three crack inclination angles are the same, as shown in Figure 6(b).A rectangular plate with an arbitrarily oriented crack. (a) Geometry and loading, (b) computational meshes.
(a)
(b)Figures7 and 8, respectively, present the normalized mode-I and mode-II DSIFs for the three considered crack angles computed by the X-IGA, sES-FEM [37], and X-FEM. Figure 6(b) is also the mesh for the X-FEM, whereas the mesh of 60 × 120 triangular elements is used in the analysis of the sES-FEM.Figure 7
Comparison of the normalized mode-I DSIFs obtained by the X-FEM, the sES-FEM, and the X-IGA for three crack inclination angles.Figure 8
Comparison of the normalized mode-II DSIFs obtained by the X-FEM, the sES-FEM, and the X-IGA for three crack inclination angles.A good agreement can be observed for the three considered crack inclination angles among the three numerical approaches.


### 5.3. Edge Crack in a Rectangular Orthotropic Plate

As a third example, a rectangular orthotropic plate with a 12 mm edge horizontal crack is considered (Figure9). The plane stress condition is considered, and the material properties of the plate are as follows: E
1
= 118.3 GPa, E
2
= 54.8 GPa, G
12
= 8.79 GPa, ν
12
= 0.083, and ρ
= 1900 kg/m3. The time step Δ
t
=
(
a
/
10
c
L
) is selected, where c
L is the wave velocity along the E
2 material-axis.Figure 9
Geometry and loads of an edge crack in a rectangular orthotropic plate.Figure10 shows the mesh with 75 × 57 elements used in the simulation. The normalized mode-I DSIFs achieved with the X-IGA (Figure 11) were compared with the reference boundary element method (BEM) (21 elements for the external boundary and eight elements for the crack) [41], the X-FEM (78 × 60 quadrilateral elements) [42], and the conventional finite element method (FEM) results (via ANSYS) [41]. Good agreement is observed among the methods. The X-IGA results more closely match the FEM and the BEM solution compared with the X-FEM results. In [42], the DSIFs are evaluated by using the domain separation integral method, which could be the reason that the X-FEM results are different.Figure 10
The X-IGA discretization for an edge crack in a rectangular orthotropic plate.Figure 11
Comparison of the normalized mode-I DSIFs among different methods.To analyze the effects of different meshes on DSIFs, five different meshes with 27 × 19, 47 × 39, 75 × 57, 97 × 81, and 121 × 101 elements are adopted to simulate the problem. The results are shown in Figure12. Accurate results can be obtained by using the X-IGA with a relatively coarse mesh.Figure 12
Comparison of the normalized mode-I DSIFs for different meshes.


### 5.4. Central Crack with Different Orientations of the Axes of Orthotropy

A rectangular orthotropic plate with a single central horizontal crack was subjected to a Heaviside step tensile distributed load, as shown in Figure13(a). The length of the crack is 2
a
= 4.8 mm and h
= 20 mm. The material properties of the plate are the same as those in the example in Section 5.3. The plane stress state is assumed, and the time step Δ
t
=
(
h
/
50
c
L
) is selected.A single central crack in a rectangular orthotropic plate. (a) Geometry and loadings, (b) computational meshes.
(a)
(b)Two different inclination angles between theE
1 material-axis and the crack-face are considered. Figures 14 and 15 show the comparison of the numerical results of BEM (24 elements for the external boundary and 10 for the crack) [41] and X-FEM (50 × 100 quadrilateral elements) [42]. General trends are in good agreement across all methods; the X-IGA results more closely match the BEM solution compared with the X-FEM; the reason for this phenomenon is the same as the above example.Figure 14
Normalized mode-I DSIFs for a single central crack in a rectangular plate.Figure 15
Normalized mode-II DSIFs for a single central crack in a rectangular plate.


### 5.5. Edge Crack in an Annular Isotropic Plate

To illustrate the capability and versatility of the X-IGA in modeling the complicated geometries, the last example deals with an annular isotropic plate under a uniform pressure. The geometry and loads of the plate are depicted in Figure16. The length of the crack is a
= 2 m; the material properties of the plate are as follows: Young’s modulus E
= 210 GPa, Poisson’s ratio ν
= 0.3, and mass density ρ
= 8000 kg/m3. The total time of the simulation is 0.004s. The plane strain state is assumed, and the time step Δ
t
= 8 × 10−5 s is selected.Figure 16
Geometry and loads of an edge crack in an annular isotropic plate.The meshes are depicted in Figure17 for 285 control points and 221 elements. To perform a comparison, this example is solved by the X-FEM, and the X-FEM meshes are shown in Figure 18 for 1292 nodes and 1221 quadrilateral elements.Figure 17
Meshes and control points. (The cross represents the control points enriched by Heaviside function; the square represents the control points enriched by crack tip enrichment functions.)Figure 18
The X-FEM meshes.The DSIFs are computed with the X-FEM and the X-IGA, and the results are shown in Figure19. It is found that the results from the X-IGA with a relatively coarse mesh are in good agreement with those from the X-FEM with a very fine mesh. The amplitude of the normalized mode-I DSIF is negative in some time ranges, which may be due to the contact between two crack-faces which is not taken into consideration in the present study.Figure 19
Comparison of the normalized mode-I DSIFs between the X-FEM and the X-IGA.


## 5.1. Semi-Infinite Crack in an Infinite Isotropic Plate

For the first example, we consider an infinite isotropic plate with a semi-infinite crack subjected to a tensile stress wave perpendicular to the crack face, as shown in Figure1. The material properties of the media are as follows: Young’s modulus E
= 210 GPa, Poisson’s ratio ν
= 0.3, and mass density ρ
= 8000 kg/m3. The tensile stress σ
0
= 500 MPa. The theoretical solution of mode-I DSIF for the stationary crack is expressed as [38]
(25)
K
I
dyn
(
t
)
=
{
0
,
t
<
t
c
2
σ
0
1
-
v
c
d
(
t
-
t
c
)
(
1
-
2
v
)
π
,
t
c
≤
t
≤
3
t
c
,
where t
c
=
(
H
/
c
d
)  (c
d: the dilatational wave speed).Figure 1
Geometry and loading of a semi-infinite crack.


### 5.1.1. Convergence Study of the DSIF versus Meshes

In this study, five meshes with 27 × 9, 47 × 17, 67 × 25, 87 × 33, and 107 × 41 elements are considered. Figure2 illustrates the typical regular mesh of 27 × 9 elements.Figure 2
Control points and elements for degree 3 X-IGA with 27 × 9 elements (the blue cross represents the control points enriched by Heaviside function; the blue square represents the control points enriched by crack tip enrichment functions; the red line represents the crack).Figure3 presents the normalized mode-I DSIFs as a function of time for the five considered meshes. In the mesh refinement process, the mode-I DSIFs obtained by the X-IGA converge well to the theoretical solution. The maximum error increases beside t
=
t
c for the coarse mesh and decreases significantly for the fine mesh. Similar observation was noted in the simulation of the singular edge-based smoothed finite element method (sES-FEM) [37]. The same simulations were conducted by other researches by using different methods. However, the same conclusion is recorded; that is, the error is lower in later stages than in earlier stages. The relative errors during 1.5
t
c
≤
t
≤
3
t
c are presented in Figure 4, in which the error is unchanged as the mesh is refined to a certain extent. In this study, the errors are the same for meshes with 67 × 25, 87 × 33, and 107 × 41 elements.Figure 3
Convergence of the normalized mode-I DSIFs versus mesh densities.Figure 4
Convergence in percentage errors for the normalized mode-I DSIFs versus mesh densities during1.5
t
c
≤
t
≤
3
t
c.


### 5.1.2. Comparison with X-FEM Results

The X-IGA results using 67 × 25 mesh elements are compared with three available X-FEM solutions, including the X-FEM using mass lumping (120 × 60 quadrilateral elements) [39], the mesh-free enriched X-FEM with crack-tip enrichments (78 × 39 quadrilateral elements) [38], and X-FEM with a new enrichment function (120 × 59 quadrilateral elements) [40]. Figure 5 shows the comparative study of the normalized mode-I DSIFs. The results obtained with X-FEM using mass lumping [39] are unstable, and the accuracy of this method is the worst recorded among the four approaches. The results obtained with the X-IGA, the X-FEM using an enrichment function, and the mesh-free enriched X-FEM with crack-tip enrichments are very accurate during 1.5
t
c
≤
t
≤
3
t
c; however, the results are relatively poor during t
≤
1.5
t
c. The X-IGA results are more accurate than other numerical results during 1.5
t
c
≤
t
≤
3
t
c.Figure 5
Comparison of the normalized mode-I DSIFs among different methods.


## 5.1.1. Convergence Study of the DSIF versus Meshes

In this study, five meshes with 27 × 9, 47 × 17, 67 × 25, 87 × 33, and 107 × 41 elements are considered. Figure2 illustrates the typical regular mesh of 27 × 9 elements.Figure 2
Control points and elements for degree 3 X-IGA with 27 × 9 elements (the blue cross represents the control points enriched by Heaviside function; the blue square represents the control points enriched by crack tip enrichment functions; the red line represents the crack).Figure3 presents the normalized mode-I DSIFs as a function of time for the five considered meshes. In the mesh refinement process, the mode-I DSIFs obtained by the X-IGA converge well to the theoretical solution. The maximum error increases beside t
=
t
c for the coarse mesh and decreases significantly for the fine mesh. Similar observation was noted in the simulation of the singular edge-based smoothed finite element method (sES-FEM) [37]. The same simulations were conducted by other researches by using different methods. However, the same conclusion is recorded; that is, the error is lower in later stages than in earlier stages. The relative errors during 1.5
t
c
≤
t
≤
3
t
c are presented in Figure 4, in which the error is unchanged as the mesh is refined to a certain extent. In this study, the errors are the same for meshes with 67 × 25, 87 × 33, and 107 × 41 elements.Figure 3
Convergence of the normalized mode-I DSIFs versus mesh densities.Figure 4
Convergence in percentage errors for the normalized mode-I DSIFs versus mesh densities during1.5
t
c
≤
t
≤
3
t
c.


## 5.1.2. Comparison with X-FEM Results

The X-IGA results using 67 × 25 mesh elements are compared with three available X-FEM solutions, including the X-FEM using mass lumping (120 × 60 quadrilateral elements) [39], the mesh-free enriched X-FEM with crack-tip enrichments (78 × 39 quadrilateral elements) [38], and X-FEM with a new enrichment function (120 × 59 quadrilateral elements) [40]. Figure 5 shows the comparative study of the normalized mode-I DSIFs. The results obtained with X-FEM using mass lumping [39] are unstable, and the accuracy of this method is the worst recorded among the four approaches. The results obtained with the X-IGA, the X-FEM using an enrichment function, and the mesh-free enriched X-FEM with crack-tip enrichments are very accurate during 1.5
t
c
≤
t
≤
3
t
c; however, the results are relatively poor during t
≤
1.5
t
c. The X-IGA results are more accurate than other numerical results during 1.5
t
c
≤
t
≤
3
t
c.Figure 5
Comparison of the normalized mode-I DSIFs among different methods.


## 5.2. Arbitrarily Oriented Central Crack in a Rectangular Isotropic Plate

We consider a rectangular isotropic plate with an arbitrarily oriented central crack subjected on the top and bottom of a uniform impact loading, as shown in Figure6(a). The length of the crack is 2
a
= 4.8 mm. The material properties of the plate are as follows: Young’s modulus E
= 200 GPa, Poisson’s ratio ν
= 0.3, and mass density ρ
= 5000 kg/m3. Three crack inclination angles are investigated, and the total time of the simulation is 20 μs. In the X-IGA, the computational mesh is independent of the crack; therefore, the computation meshes for the three crack inclination angles are the same, as shown in Figure 6(b).A rectangular plate with an arbitrarily oriented crack. (a) Geometry and loading, (b) computational meshes.
(a)
(b)Figures7 and 8, respectively, present the normalized mode-I and mode-II DSIFs for the three considered crack angles computed by the X-IGA, sES-FEM [37], and X-FEM. Figure 6(b) is also the mesh for the X-FEM, whereas the mesh of 60 × 120 triangular elements is used in the analysis of the sES-FEM.Figure 7
Comparison of the normalized mode-I DSIFs obtained by the X-FEM, the sES-FEM, and the X-IGA for three crack inclination angles.Figure 8
Comparison of the normalized mode-II DSIFs obtained by the X-FEM, the sES-FEM, and the X-IGA for three crack inclination angles.A good agreement can be observed for the three considered crack inclination angles among the three numerical approaches.


## 5.3. Edge Crack in a Rectangular Orthotropic Plate

As a third example, a rectangular orthotropic plate with a 12 mm edge horizontal crack is considered (Figure9). The plane stress condition is considered, and the material properties of the plate are as follows: E
1
= 118.3 GPa, E
2
= 54.8 GPa, G
12
= 8.79 GPa, ν
12
= 0.083, and ρ
= 1900 kg/m3. The time step Δ
t
=
(
a
/
10
c
L
) is selected, where c
L is the wave velocity along the E
2 material-axis.Figure 9
Geometry and loads of an edge crack in a rectangular orthotropic plate.Figure10 shows the mesh with 75 × 57 elements used in the simulation. The normalized mode-I DSIFs achieved with the X-IGA (Figure 11) were compared with the reference boundary element method (BEM) (21 elements for the external boundary and eight elements for the crack) [41], the X-FEM (78 × 60 quadrilateral elements) [42], and the conventional finite element method (FEM) results (via ANSYS) [41]. Good agreement is observed among the methods. The X-IGA results more closely match the FEM and the BEM solution compared with the X-FEM results. In [42], the DSIFs are evaluated by using the domain separation integral method, which could be the reason that the X-FEM results are different.Figure 10
The X-IGA discretization for an edge crack in a rectangular orthotropic plate.Figure 11
Comparison of the normalized mode-I DSIFs among different methods.To analyze the effects of different meshes on DSIFs, five different meshes with 27 × 19, 47 × 39, 75 × 57, 97 × 81, and 121 × 101 elements are adopted to simulate the problem. The results are shown in Figure12. Accurate results can be obtained by using the X-IGA with a relatively coarse mesh.Figure 12
Comparison of the normalized mode-I DSIFs for different meshes.


## 5.4. Central Crack with Different Orientations of the Axes of Orthotropy

A rectangular orthotropic plate with a single central horizontal crack was subjected to a Heaviside step tensile distributed load, as shown in Figure13(a). The length of the crack is 2
a
= 4.8 mm and h
= 20 mm. The material properties of the plate are the same as those in the example in Section 5.3. The plane stress state is assumed, and the time step Δ
t
=
(
h
/
50
c
L
) is selected.A single central crack in a rectangular orthotropic plate. (a) Geometry and loadings, (b) computational meshes.
(a)
(b)Two different inclination angles between theE
1 material-axis and the crack-face are considered. Figures 14 and 15 show the comparison of the numerical results of BEM (24 elements for the external boundary and 10 for the crack) [41] and X-FEM (50 × 100 quadrilateral elements) [42]. General trends are in good agreement across all methods; the X-IGA results more closely match the BEM solution compared with the X-FEM; the reason for this phenomenon is the same as the above example.Figure 14
Normalized mode-I DSIFs for a single central crack in a rectangular plate.Figure 15
Normalized mode-II DSIFs for a single central crack in a rectangular plate.


## 5.5. Edge Crack in an Annular Isotropic Plate

To illustrate the capability and versatility of the X-IGA in modeling the complicated geometries, the last example deals with an annular isotropic plate under a uniform pressure. The geometry and loads of the plate are depicted in Figure16. The length of the crack is a
= 2 m; the material properties of the plate are as follows: Young’s modulus E
= 210 GPa, Poisson’s ratio ν
= 0.3, and mass density ρ
= 8000 kg/m3. The total time of the simulation is 0.004s. The plane strain state is assumed, and the time step Δ
t
= 8 × 10−5 s is selected.Figure 16
Geometry and loads of an edge crack in an annular isotropic plate.The meshes are depicted in Figure17 for 285 control points and 221 elements. To perform a comparison, this example is solved by the X-FEM, and the X-FEM meshes are shown in Figure 18 for 1292 nodes and 1221 quadrilateral elements.Figure 17
Meshes and control points. (The cross represents the control points enriched by Heaviside function; the square represents the control points enriched by crack tip enrichment functions.)Figure 18
The X-FEM meshes.The DSIFs are computed with the X-FEM and the X-IGA, and the results are shown in Figure19. It is found that the results from the X-IGA with a relatively coarse mesh are in good agreement with those from the X-FEM with a very fine mesh. The amplitude of the normalized mode-I DSIF is negative in some time ranges, which may be due to the contact between two crack-faces which is not taken into consideration in the present study.Figure 19
Comparison of the normalized mode-I DSIFs between the X-FEM and the X-IGA.


## 6. Conclusions and Prospects

The NURBS-based extended isogeometric analysis (X-IGA) was extended to investigate the dynamic fracture behavior of stationary cracks in isotropic/orthotropic media. A corresponding dynamic X-IGA model was developed, the Newmark time integration scheme was used to achieve the dynamic response, and the DSIFs were evaluated with the contour interaction integral technique. Numerical results indicate that accurate results can be obtained by using the X-IGA with a relatively coarse mesh, and the X-IGA is suitable to model the cracked complicated geometries.The NURBS-based isogeometric analysis has several drawbacks including the handling trimmed geometries and the local refinement. The isogeometric analysis based on T-splines can effectively overcome these drawbacks [43]. A key advantage of the extended isogeometric analysis over isogeometric analysis is the ability to model crack propagation. Therefore, the application of extended isogeometric analysis based on T-splines to dynamic crack propagation analysis in isotropic/orthotropic media is an area of considerable promise. On the other hand, the smoothed finite element method (SFEM) developed by combining the smoothing technique with the finite element method is an efficient and accurate numerical simulation tool for the dynamic fracture problems [44]; therefore combining the smoothing technique with the X-IGA for fracture analysis is another area of considerable promise.



---

*Source: 725795-2014-04-08.xml*