User login


You are here

Crack Propagation

atmaca's picture


I have an investigation on Crack Propagation.

How can i predict the path of a crack.

Please help me!

atmaca's picture

I have a FORTRAN program and i can calculate u and v (x-y displacements) and stress values using this program.

I wonder which parameter(among displacemnt or stress values) i have to use for the prediction of crack path.

Or is it necessary to calculate K values



Roberto Ballarini's picture

You have received many replies to your question. Perhaps the quickest way for you to do what you want to do is download the program Franc-2d from the Cornell Fracture Group website.

 This program allows you to grow cracks according to the criteria mentioned in previous replies; it automatically remeshes propagating cracks, and provides accurate solutions.

 Good luck.

Dear Atmaca;

I you want, please share that code  to

We will discuss about that, my research is same you, I must modify the  code C++ for predicting the crack propagation under dynamic loading.

Very happy for your cooperation, Goodluck 

Zhigang Suo's picture

In a brittle, homogeneous, isotropic material, crack propagates along a mode I path. That is, KII = 0 at the tip of the crack as the crack advances. This path-selection criterion seems to be consistent with a number of experimental observations. One example that comes to mind is


To predict crack path by numerical calculation, one has to calculate the crack tip field accurately, and then use the mode-I criterion or some other criteria to advance the crack by a small amount. One then repeats this process incrementally.

The essential requirements for the numerical methods are then

  1. Resolve stress field around the crack tip.
  2. Track the crack path.

These requirements will challenge traditional finite element method. I have seen people use boundary element methods and meshfree methods. Perhaps experts of these methods can give you the references. Here is a paper in which I was involved:

J. Liang, R. Huang, J.H. Prévost, and Z. Suo, Evolving crack patterns in thin films with the extended finite element method. International Journal of Solids and Structures. 40, 2343-2354 (2003).

This paper applies the meshfree method, and gives some of the references.

Edit:  meshfree method above should be replaced by XFEM.  See the comment by YingLee

Ying Li's picture

You can see the follow URL:

The XFEM is pretty for the Crack Propagation calculation.

allspot's picture

Good share bro thanks ....

N. Sukumar's picture


A few more details; Zhigang's response is concise and clear.  Brittle crack propagation criteria come in different flavors but they all lead to similar end  results, whether one uses the maximum hoop stress, maximum energy release rate or the maximum strain energy density criteria. In 2D, to apply any of these, one needs to determine K_I and K_II from the numerical solution (displacements, strains, and stresses are first computed in the domain).  For FEM and method of that ilk, the domain form (Li, Shih, Needleman; Shih, Moran and Nakamura; Nikishkov and Atluri; Yau, Wang and Corten, etc.) of the interaction integral is well-suited to extract the mixed-mode SIFs. For pure mode I, the above reduces  to the domain form of Rice's contour J-integral. The above are spelled out in many of the x-fem papers and also in the FEM literature: derivation of the domain form of the  interaction integral is provided in Moes, Dolbow and Belytschko and some  of the above details/references appear in the x-fem paper with Prevost.

In 3D, the problem is a lot more demanding. One of the outstanding issues from a computational perspective is how best to extrack mixed-mode SIFs/energy release rate along a crack front in 3D for Galerkin methods.  In 3D, are there better alternatives than domain integral formulations for SIF computations (accurate and reasonably fast)?

Henry Tan's picture

For an even harder and interesting problem, how to simulate a 3D body with lots of interacting cracks?

N. Sukumar's picture

Yes,  that's a tough one; just dealing with 2 planar cracks posed many challenges and (numerical) issues.  In all cases, need pretty good/accurate SIF calculations on the crack front so that the front remains smooth and is less prone to wreaking havoc in the subsequent `time' steps (quasi-static).

The boundary integral equation/boundary element method (BIE/BEM), especially the one accelerated with the fast multipole method (FMM), is pretty good in modeling multiple crack interactions and growths in both 2-D and 3-D. In the BEM, only the boundary (e.g., crack surfaces) need to be discretized and more accuracte values of SIFs can be obtained. 2-D models with thousands of cracks and 3-D models with hundreds of cracks have been solved successfully with the fast multipole BEM. See two recent papers:

  1. Wang, P. and Z. Yao (2006). "Fast multipole DBEM analysis of fatigue crack growth." Computational Mechanics 38: 223-233.
  2. Yoshida, K., N. Nishimura and S. Kobayashi (2001). "Application of fast multipole Galerkin boundary integral equation method to crack problems in 3D." International Journal for Numerical Methods in Engineering 50: 525-547.

Julian J. Rimoli's picture

A paper by Ortiz & Pandolfi proposes 'Finite-deformation irreversible cohesive elements for three-dimensional crack-propagation analysis'

The paper can be found at:

I think the method proposed is very powerful for simulating 3D crack propagation problems and its integration into standard FEM codes is quite straightforward (you just have to implement the tangent stiffness matrix).

I hope it helps.



atmaca's picture

Thank you for all repplies


Roberto Ballarini's picture

You have received many replies to your question. Perhaps the quickest way for you to do what you want to do is download the program Franc-2d from the Cornell Fracture Group website.

 This program allows you to grow cracks according to the criteria mentioned in previous replies; it automatically remeshes propagating cracks, and provides accurate solutions.

 Good luck.

Dear Nihat,

I would be pleased to help you with your problem, please contact me on bordas AT civil DOT gla DOT ac DOT uk I will give you a computer code based on XFEM which deals with crack propagation.  

Dr Stephane Bordas

You can also take a look at this paper in IJNME, which should help you understand the code.  

An extended finite element library (p n/a)
Stéphane Bordas, Phu Vinh Nguyen, Cyrille Dunant, Amor Guidoum, Hung Nguyen-Dang
Published Online: 5 Jan 2007
DOI: 10.1002/nme.1966

Dr Stephane Bordas

atmaca's picture

Dear Stephane Bordas

I have send e-mail at your gmail adress.

Thank you.


Mahdi Kazemzadeh's picture

Dear Nihat,

This free programme will help you to do mesh configuration that you have in your mind. Actually first of all you define initial cell vectors which are about your sample under crack, by giving the number of steps and time lattice you will find this programme changing your initial position by time and the other factors which you have and you can improve. This file attached below will be additional aid to define initial position of your crack or even use MOLDY to do it for you, please feel free to contact me if you will need more information.

Additional Information about cell vectors

Good Luck. 

Mahdi Kazemzadeh

HI MR Stephane

please do you mind if you send me the code dealing with the crack?, i have the same problem

as Mr Nihat.


Can you please contact me by email? Thanks,   

Dr Stephane Bordas



May you email me a copy of your XFEM code? 

I want to see if I can use it to do simulations of some fluid shock wave problem.

Thank you!




For mode I  and mixed-mode I/II elastic-plastic problems, you may try CTOD criterion. The following is a recent paper by my advisors and me:

 Lan, W., Deng, X., Sutton, M.A., Three-dimensional finite element simulations of mixed-mode stable tearing crack growth experiments, Engineering Fracture Mechanics (2007), doi: 016/j.engfracmech.2006.12.026.

Mahdi Kazemzadeh's picture

Dear Nihat,

You have got some solutions which they all are helpful. I had the same issue before and here is some points which I think will be useful later and if you want to continue working on this issue,

1) There is free MD (molecular dynamics) simulation programme on web which you can download from link below using FTP server. The name is MOLDY and it is very simple to work with. In order to solve your problem you can use just some sub. of this programme after compiling it. That is written in C and later you can develope it for your further multi scale modellings if you want. It is possible to use parallel computation systems for this purpose. If you couldnt find the place to develope this command lines for your issue, feel free to email me.

 A GENERAL-PURPOSE MD CODE  or type Moldy in Google and click on the name above.

2) you can write system specification file for your initial conditions that exist in your issue and then by running this model you will get some results.

3) Defining cell vectors, temperature, moleculs diameter and your modelling steps will enable you to simulate propagation and solve your problem even numerically. It is easy to work with.

4) Considering stress distribution nearby crack tip, for brittle isotropic material you have to come up with this (K). It varies with the intensity of stress around your tip and your mode which can be opening mode (I). It depends on your crack shape and physics of your material also. MOLDY will help by a right assumptions which you will do.

5) Some scientist like, Prof. gumbsch, Prof. Gao, Dr. Buehler even have worked on hyperelasticity which occures in crack tips. Here is the list of papers helped me come up with this issue somehow:

- Molecualr dynamics investigating of dynamic crack stability, p. Gumbsch, 1997.

-A molecualr dynamics investigation of rapid fracture mechanics, F. Abraham, 1996

-Atomistic Aspects of Crack Propagation in Brittle Materials, Elefterios Lidorikis and ..., 2002

-Dynamical fracture instabilities due to local hyperelasticity at crack tips, Dr. Buehler, 2006

They are not cited in a normal way!! if you couldnt find them I will send you the link from Scopus. By devoting more time in this stage of your work, modelling will more enjoyable later. Please feel free to contact me if I can help anyway.

Good Luck,

Mahdi Kazemzadeh

Dear All,

The direction of the crack growth depends on the current DSIFs and the crack velocity or only denpends on the crrent DSIFs.

Because i have read the method to predict the direction of crack growth in the book of Broek (Elementary engineering Fracture Mechanics, page 375- 376) and i saw that the direction of crack growth only depends the current DSIFs.

But in the paper (The time-domain DBEM for rapidle growing cracks of P. FEDELINSKI* AND M. H. ALIABADIR,  INTERNATIONAL JOURNAL FOR NUMERICAL METHODS IN ENGINEERING, VOL. 40, 1555-1572 (1997)), i have seen that the  direction of crack growth depends on DSIFs and velocity of crack tip.

So, please give me comment about that problem.



Hi, I also have written a FORTRAN 90 crack propagation FEM code but for 2D elastic case. I use displacement exptrapolation technique to calculate the stress intensity factor KI and KII. Currently I am trying to calculate J-integral using equivalent domain integral but somehow need to create mesh domain (ring) around the crack tip first. I would be happy if somebody can help me out.

If you mean that you want to produce a closed node-path around the crack tip which is idepedent of the radius-distance in a local  CS from the crack tip, (you should rettrive -about-the same results ) I have a code in ANSYS APDL for doing this.It produces the path and then calculates the J integral using the values estimated in the nodes of the path.

From my reading calculating J-integral using close path line integral has disadvantage as SIF for mode I and II for example cannot be separated unlike using EDI.

Octavio A. González Estrada's picture

Hi, if you're looking to obtain both SIFs it's better that you use an Interaction Integral. As with EDI, you'll need to define a plateau function 'q' that takes values of 1 for nodes within a circle of user-defined radius centered on the crack tip, and 0 for the remaining nodes. Defined this way, you evaluate the integral on elements around the crack tip where the q-function has gradient different from zero, i.e. the ring you mention.


I this the best way is you read some code which related to crack propagation. They are free on Internet. If you want to get the nodes around crack tip to use for calculate the dynamic J-integral, you should choose a circle which its center is crack tip point and on this circle we will some nodes, these nodes will be uesed to calculate the integration under line integral mark of J-path formular. And the other nodes inside that circle will be used to calculate the integration under area integral mark of J-pth formular. When you change the radius of the circle, those datum in above procedure will change but its value is not change. This is called the independent J-path.

Good luck

Fracture man...


i have the same probleme like Mr Nihat

please could you give me more details !


allspot's picture

Dear, Mr Fracture Man Cool

where I can get that source code Sir

Maybe you can share .... 

I need  it to validation result output from ansys ...


Thanks before and greetings super 

dear atmaca

you need to things to make crack propagatio

1-criteria for the crack propagation direction

     i know some like


a) The Maximum Energy Release Rate Criterion

b) Maximum Circumferential (Tangential) Stress Criterion (MTS-Criteria)

c) Minimum Strain Energy Density Criterion (S-Criterion)

(d) The J-criterion

e) Maximum Dilatational Strain Energy Density Criterion (T-criterion)


f) Maximum Dilatational Strain Energy Density Criterion (Y-Criterion)

g) Vector Crack Tip Displacement (CTD) Criterion

h) Tangential Stress Factor and Tangential Strain Factor Criteria


i) Maximum Tangential Strain Criterion

k) Ip-Criterion


Among the aforementioned criteria, the maximum tangential stress and the

minimum strain energy density criteria are widely used in mixed mode crack

growth studies.

 2-u need  a low  to calculate the number of cycles for each crack icr

i make my models by ANSYS i insert code for remesking in it and for fracture parameters calculation


Im a student of rock mechanic engineer(master)


My thesis is about crack propagation in rock slope (brittle materials) with mathematica 6Before use G ore S ore sikma  critration I need a failoure critration for start of Crack appearance and this is a problem.I need a sample of mathematica based on crack propagation.



   Do you have some doc regarding those criteria ?

Tks !



   i m a 2nd year student doing project on fracture mechanics. i want to find the SIF for mutiple cracks in same model. can anybody tell me how to do that? 

damoonmoetamedi's picture

Hi all,

I'm working on crack propagation through XFEM. I chose  Maximum Circumferential Stress Criterion and in formulation, we need to know that KII is positive or negative. (θc=2*arctan(0.25(KI/KII-sign(KII)*sqrt((KI/KII)^2+8)))).

(Because I used seperation J-integral method) 

1. How could I find that the KII is positive or negative?

2. Is there any other criterion that you recommend?


N. Sukumar's picture

You can use the expression for the crack growth direction that was
suggested by Suo and appears on page 7531 in S and Prevost, IJSS, V40,
7513-7537, 2003. Only the ratio KII/KI appears; do not need to know the
sign of KII.


damoonmoetamedi's picture

Dear Dr.Sukumar,

Thanks for hints, but I think if we use the formulation, the critical direction always will be negative (considering resullt of arctan between ∏/2 & -∏/2).

Because this angle is calculated to local direction of crack, it will be sum up with previous crack angle and so, direction of crack won't be correct, unless similar to Maximum Circumferential Stress criterion, KII shows positive and negative values.

And, What is the meaning of negative KII? This fracture mode shows sliding of crack surfaces, how could we interpret this problem?



Rui Huang's picture


No matter what criterion you use, the direction of crack growth would depend on the sign of KII. Once the coordinate is set, the sign of KII depends on the sign of the shear stress ahead of the crack tip, which can be either positive or negative. Note that KI is always positive. For a straight crack, it kinks upward when KII is negative and kinks downward when KII is positive, in a standard crack coordinate. In XFEM, one may calculate KI and KII by using the interaction J-integral method as outlined in the Appendix of this paper:

J. Liang, R. Huang, J.H. Prévost, and Z. Suo, "Evolving crack patterns in thin films with the extended finite element method".
International Journal of Solids and Structures
,40, 2343-2354 (2003).

Hope it helps.


damoonmoetamedi's picture

Dear Huang,

Your hint seems to be really helpful. I'm going to implement it in my source and I think the result will confirm your idea.

Today, I have used Maximum energy release rate by mean of J= J1cosθ + J2sinθ. (θc = arctan(J2/J1)).  The result was in agreement with other available examples.

If you have any idea about this method, I will be glad to hear from you.

Great thanks,


Rui Huang's picture

I have no experience using the maximum energy release rate criterion in FEM or XFEM, but am aware of some analytical approaches using this criterion. As far as I know, it requires calculation of energy release rate of a virtual extension of the current crack in various directions, which may not be convenient for numerical implementation. As for your equation, I am not sure how you define J1 and J2. I am not aware of a solution in this form. In any case, the angle θc can be either positive or negative, and thus one of your J1 and J2 must be allowed to take negative values.



damoonmoetamedi's picture

Dear Huang,

I am using following formulationto calculate KI and KII:


Also, for crack speed I use (( Vcr = (1-Greleased/Gcritical)*Cs)). 

Now, when crack propagate with Vcr, the J1 decreases but it gets bigger values rather than previous steps and sometimes negative value (J2 too), which I dont think it is logical.

Please give me some advice. 


Are you trying to say that the J is a vector?


Die ersten hundert Jahre sind erfahrungsgemaess die Schwersten!
Alles Gute!

Rui Huang's picture


Looks like you are dealing with a dynamic crack growth problem, with which I am less familiar. Other imechanicians, like K. Ravi-Chandar and Young Huang , may have better ideas about the dynamic problem. 


damoonmoetamedi's picture

Dear Huang,

Thanks for your care and responds.

 Best wishes,


as's picture

Maybe this reference can be of interest:  A. Salvadori: “A plasticity framework for (linear elastic) fracture mechanics” - Journal of the Mechanics and Physics of Solids, 56 (2008) 2092–2116. Hope it helps, Alberto






TO LEARN ABAQUS, come this web:

Dear alls,

I have some points which i not yet understood, can you give me some comments:

1. Maybe you know the Book of Anderson (Fracture Mechanics), 3rd edition. In the chapter 4, page 183, example 4.1, i want to know why he use the longitudinal wave speed Cl, because you know, the dilatational wave speed Cd bigger longitudinal wave speed Cl, right?

So dilatational wave in the material body will come to crack tip sooner Cl.

2. Did you know the paper or the other source say about the sinificant effect of reflected stress wave on the local crack-tip fileds?

Please help me!!

Thanks and Regards,,


If you think of the crack tip as a point dislocation source, then there are essentially two types of elastic waves that propagate.  These are the P-wave (pressure wave or primary wave) that shows up first in an accelerometer trace and the S-wave (shear wave or secondary wave) which shows up later.

The P-wave is also called the "longitudinal wave" because the direction of the particle motion is the same as the direction of propagation.  The velocity of this wave is  c_p = sqrt{(lambda+2 mu)/rho}.  This wave is also called the "dilational" wave.  See for example "Waves and Imaging Through Complex Media ". 

A good book for the basics is "Wave Propagation in Elastic Solids" by Jan D. Achenbach.  Another that is good for applied wave propagation is "Quantitative Seismology - Vol 1" by Keiti Aki and Paul Richards.

-- Biswajit 


 I am no engineer, but live in Japan and see many kinds of cracks every day, and sometimes new ones after an earthquake. I am intrigued by how cracks appear to radiate out from the corners of squares of concrete on the ground, or from the corners of window and door frames. I guess predicting where cracks go could largely depend on being able to predict where they start.

It gets more scary when we see news reports on the abundance of cracks in road tunnels, train tunnels, and the like. It is not clear what keeps all the infrastructure together, apart from hope.


The Research Cooperative
(NPO) - an online meeting place for researchers,science writers, editors, translators, publishers and others

Dear Peter,


1. Roughly speaking, cracks form in the regions of stress concentration. Even after they are formed, they would propagate in the region of high local stresses or strain energy. Some of the expert discussions are about deciding which measures of local stresses/energy could best characterize the aforementioned "concentrations," in different contexts. (At other times, the crack-tip may be taken to be infinitely sharp in which case the stress field is singular, and so, the criteria of crack-growth to be applied have to be global in nature.)

2. The wall-cracks seen radiating from the corner of a window or a door are not necessarily deep. They may exist only in the outer plaster coating. If so, they pose no threat to structural integrity.

And, BTW, these cracks begin radiating from the corners because corners are the places of stress concentrations. In fact, any sudden change of geometry would act as a stress concentrator.

This form of cracking (in outer thin layers), is actually used as a means of experimental measurement of stresses. The technique is called "brittle coatings." Browse through the books on experimental stress analysis, or do a search on the Internet. For example, see the following URLs: (It gives a nice photograph) (These class notes explain the thing a bit.)

3. As to the cracks you see on tunnels, I am guessing that they too must have formed in the outer plaster coating only. However, one would need to have a look at a few photographs or so before one could tell.... But I am sure in Japan there are concerned engineers anyway.

4. Overall, please note. Stress fields, being tensor fields, are mathematically more complicated than many vector fields such as, e.g., EM fields, or velocity fields in the laminar flow (as traced out by, say, ink streaks in clear water) are.

freedom teaching and research


Im a student of kntu aerospace in iran(master)

My thesis is about crack propagation  with using xfem code in fortran

 I need a help for that code.


Yancheng ZHANG's picture

hello, every one

recently,I need to develop model on the material of Ti-6Al-4v with Abaqus, I will use both Johnson-cook and shear failure model,
 i am confused on the parameters for shear failure:
ks, i didn't find this value on internet;
the fracture energy for damage evolution, the value should be the same with Johnson-cook,help example madkes the same value, but i think it should be KII value for shear fracture,not KI.

  I don't know if I have expressed clearly,  could you help me ?

recently, i need to do a project it title is "Application of fracture mechanics to predict failure of a cracked structure".

can anyone tell me how 2 predict failure of cracked structure, and what parameters i need to consider to predict crack.......


Dear All,

How we can do Crack modeling of RC beam in Ansys?

How  we can investigate Crack-path of RC beam in Ansys ?

Pls send me if you have pdf file( step by step learning file)!

Pls help me!

my email is




Hi everyone

I am new user in Abaqus, i read most of the posts related to the modelling the crack but unfortunately i didnt get a clear procedure, maybe it is due to my weakness in this software.

i want to model a crack in abaqus and analysis the crack propagation, actually i am working on piezoelectric material, but i want to know,  can you tell me different steps to model a crack and then how to analyze the propagation on it.

Please help me is this case.  

Hallo everyone,

 i m using the Subroutine UDMGINI in Abaqus, i m looking for  a method to get the crack-tip coordinate in Subroutine in Fortran.

is it possible?

i will be really grateful if anybody can help me!


Hello everyone,

I  want to perform o droptest of human Skull having crackin it.

Can I perform the simulation in Abaqus/explicit?

As far as I came to know Xfem in abaqus cant be carry out with dynamic loading.

Reply asap 


Subscribe to Comments for "Crack Propagation"

Recent comments

More comments


Subscribe to Syndicate