# calculation of flexural stress in 2D condition

Hi,

I am trying to calculate the bending stress in four-point bending test in 2D condition, but some concepts make me confused. Theoretically, the stress is the tensile stress of the bottom fiber in the beam. I know that in 3D condition, the flexural stress can be easily obtain by FL/bd^2, where b is the width of the beam cross-section. But since in 2D condition, there is no thickness for the beam or the thickness is the diameter of the sphere particle? how do I get the stress based on the applied force F?

I am thinking that instead of get stress by F, I can directly obtain the stress by using bodyStressTenso

Regards,

Steve

## Question information

- Language:
- English Edit question

- Status:
- Answered

- For:
- Yade Edit question

- Assignee:
- No assignee Edit question

- Last query:
- 2021-03-15

- Last reply:
- 2021-03-15

Karol Brzezinski (kbrzezinski) said : | #1 |

Hi Steve,

Please note that 2D solutions are used either for plane streets or plain strain. In such a case the load is actually 'linear pressure' and the force is distributed over assumed thickness (e.g. unitary). Thanks to this, you can you use your formula, by assuming the thickness is known. In may opinion there is no good answer for your question. I would suggest to take step back and rethink whether 2D is a good idea.

If you want to stick to 2D you could also 'measure' forces in the symmetry o plane and 'smear' them over a known cross-section of the sample that you are modeling. By smearing I mean you could calculate the bending moment and apply it to the theoretical beam.

Best wishes,

Karol

steve (steve-cash) said : | #2 |

Thanks Karol,

Do you mean that the force (f) I obtained from 2D is already a linear load and so that the total force over the assumed thickness (b) would be f*b, which makes the formula become sigma=f*

And I don't really understand the second part like what is the symmetry o plane, could you give more elaboration on this?

Best regards,

Steve

Karol Brzezinski (kbrzezinski) said : | #3 |

Hi,

sorry, forget the second part. It would be good for the statically indeterminate beam. Since you have a four-point bending test, you know the moment already.

Regarding the first part, I meant that in actual 2D formulation you do not have to worry about thickness.

The load that you obtain (apply?) in the simulation has a force dimension, thus I think it is related to some thickness. Whatever thickness you assume, you should treat the force as integral over this thickness. For example, if your sample has b = 40 mm, and your model has b' = 1 mm, you should apply only f'=f*b'/b = f/40.

In general, b' doesn't have to correspond to the thickness of your sample, you could also assume b'=b. This would mean that your sample is 'compressed' to the thickness of the model. In such a case, however, the stiffness considerations become harder.

In my opinion, assuming thickness to solve this problem is quite intuitive, and could be a good solution in some cases. However, it is probably oversimplified. It would be useful to know the purpose of the simulation and what is the advantage of using DEM? Is this advantage affected by the 2D assumption?

Besides all that 'fundamental questions', I think that both solutions proposed by you are quite decent (if you stick to 2D simulation ;) ). Only, remember to apply the force corresponding to the assumed thickness.

Cheers,

Karol

steve (steve-cash) said : | #4 |

Thanks Karol,

I am sorry if I my statement was confusing, let me just explain it again. So I've already obtained the flexural stress & strain form the experiments with a known dimension b=100mm. And then I construct a 2D model in yade in order to simulate this experiments and see the crack pattern. However, during the calibration process, I have some problem in matching the 2D simulation results to the experiment data. Actually, I am not applying force in both experiment and simulation, but control the loading rate of the loading plate and obtain the force that acting on the plates. So in this case, is the force obtained in the simulation supposed to multiply a scale 100/b' (100 is the thickness in the experiment, b' is the thickness of the model) to restore the force to the 3D condition and then use the formula to calculate the stress? If so, then what would be the thickness of the model? is it the average diameter of the particle? or is there no thickness at all because it is 2D simulation and both force and displacement in the out-of -plane direction has been restricted?

And I test both methods in simulation which are by formula sigma=fL/d^2 and by directly get from particles:

###

sigma=bodyStres

tensile_

###

Although the shape of curves are similar, but the peak flexural stress obtained directly from particles are over two times larger than that of formula. I am really confused where this difference come from if both methods are good enough.

Best regards,

Steve

Jan Stránský (honzik) said : | #5 |

> both force and displacement in the out-of -plane direction has been restricted?

This is not possible. In all mechanics, you know/prescribe/

In plane stress, you know the out-of-plane stress is zero, but you cannot say anything about out-of-plane strain.

In plane strain, you know the out-of-plane strain is zero, but you cannot say anything about out-of-plane stress.

....

> If so, then what would be the thickness of the model? is it the average diameter of the particle?

sounds reasonable

> or is there no thickness at all because it is 2D simulation

this is another approach. Then your model is "per unit thickness".

So for transfer model -> experiment, multiply your model data y real beam thickness.

Or for experiment -> model, divide real data by beam thickness

> the peak flexural stress obtained directly from particles are over two times larger than that of formula. I am really confused where this difference come from if both methods are good enough.

Have you done the "porosity correction" for bodyStressTensors [1]?

cheers

Jan

[1] https:/

steve (steve-cash) said : | #6 |

Hi, Jan

Thanks for the correction. I am a bit unclear why there can be other approaches for the transfer. I mean the thickness of the model must be either 0 or some specific value(eg. the average diameter). So the unit of the force I obtained from the simulation is Newton per unit thickness or just newton?

And yes I forgot the porosity correction. But after adding the porosity which is 1-(area of the all particles/ sample area) in 2D, and multiply stress by the 1-porosity, the value is still 1.7 times larger than from formula since the porosity is only around 0.2. Please correct me if there's something wrong.

Best regards,

Steve

Jan Stránský (honzik) said : | #7 |

> I am a bit unclear why there can be other approaches for the transfer

it is a model. For the same real situation, you can have many different models...

> I mean the thickness of the model must be either 0 or some specific value(eg. the average diameter). So the unit of the force I obtained from the simulation is Newton per unit thickness or just newton?

Yes.

For real thickness, you get real forces, Newtons.

For "zero" thickness (I would not call the thickness zero, but rather something like planar or 2D case, since the thickness is not zero, it is just not important for the solution itself), you get force per unit thickness, Newton per unit thickness.

> But after adding the porosity which is 1-(area of the all particles/ sample area) in 2D

> the porosity is only around 0.2

"area" porosity, or "volume" porosity?

As I mentioned, bodyStressTensors are for 3D, are used (and assumes working with) with volumes, not areas.

If you have significant differences, in particle sizes, "area" correction is wrong, since the particle "volume" and "area" are not proportional.

Maybe you can use the same math theory, only adjusted to the 2D case.

For (almost) static case, the "stress" approach should fit well with the theory.

Anyway, I think we came to the place, where a MWE [2] is needed for reasonable continuation.

If you prepare a code, please focus on the "M=minimal" part (both from code size and running time point of view)

cheers

Jan

steve (steve-cash) said : | #8 |

Hi Jan,

Sorry that I thought the simulation was 2D so I used the 2D porosity. But I forgot that stress is also obtained in 3D.

1. After correction, the peak values do become closer to each other, but there is still some gap as well as the shape.

2. the error in volume calculation is that I didn't use aabbExtrema() to calculate the volume because if I add cement particle before the roller(moving line78,79 to line 69) , an error called 'Assertion `!swap' failed' occurred.

3. I tried the TesslationWrap() but the value is way much smaller, maybe it's due to the reason that the function also considers the rollers and supports?

line 1-59 is just making 2D packing, line 60-156 is 4-point bending model

Many thanks,

Steve

#######

## compact 2D packing

from yade import pack,plot

nRead=readParam

num_spheres=2000,

compFricDegree = 30,

key='_

unknownOk=True

)

from yade.params import table

num_spheres=

key=table.key

targetPorosity = 0.55

compFricDegree = table.compFricD

finalFricDegree = 30

rate=-0.02

damp=0.2

stabilityThresh

young=5e6

mn,mx=Vector3(

O.switchScene(

O.materials.

O.materials.

walls=aabbWalls

wallIds=

sp=pack.

sp.makeCloud(

O.bodies.

triax=TriaxialS

maxMultiplier=

finalMaxMultip

thickness = 0,

stressMask = 7,

internalCompac

)

newton=

O.engines=[

ForceResetter(),

InsertionSortC

InteractionLoop(

[Ig2_

[Ip2_

[Law2_

),

GlobalStiffnes

triax,

TriaxialStateR

newton

]

triax.goal1=

triax.goal2=0

while 1:

O.run(1000, True)

unb=unbalance

print ('unbalanced force:',unb,' mean stress: ',triax.meanStress)

if unb<stabilityTh

break

sp=SpherePack(); sp.fromSimulation()

O.switchScene()

## constuct 4-point bending model

import numpy as np

from yade import pack,plot

## define parameters

strainRate=-.02

sphereRadius=.0005

young=2.8e9

sigmaT=1e6

frictionAngle=

relDuctility=17

omegaThreshold=

epsCrackOnset=.004

damping=.6

damLaw = 0

relFuzz=0.3

factor=1.5

isoPrestress=0

poisson= 0.1

pred=pack.

idConcrete=

#

idRoller=

## create top roller clumps

TL=O.bodies.

TR=O.bodies.

## creat bottom roller

O.bodies.

O.bodies.

#

spp=filterSpher

cement=

#### calculate the porosity for volume

volume=0

for i in cement:

volume+

porosity=

print(porosity)

T=[O.bodies[i] for i in cement if (O.bodies[

initial_dis= T[0].state.

####grab the particle at the bottom fibre

p=[]

for i in cement:

if (O.bodies[

p.append(i)

O.engines=[

ForceResetter(),

InsertionSort

Bo1_

]),

InteractionLoop(

[

Ig2_

],

[

Ip2_

Ip2_

Ip2_

],

[

Law2_

Law2_

]

),

GlobalStiffne

NewtonIntegra

CpmStateUpdat

PyRunner(

PyRunner(

PyRunner(

]

O.bodies[

O.bodies[

O.step()

bo1s.aabbEnlarg

ig2ss.interacti

def initiateStrain():

O.bodies[

O.bodies[

def finish():

if (O.bodies[

O.pause()

pass

def addPlotData():

f_up=

### directly obtain stress on the bottom fibre

TW=Tesselatio

TW.setState()

TW.computeVol

sigma=

tensile_

tensile_stress2= sum(sigma[

plot.

plot.plots=

plot.plot()

O.run()

## Can you help with this problem?

Provide an answer of your own, or ask steve for more information if necessary.