Page 2 of 2

Re: Stress Solver - Pre-strain is not applied?

Posted: 04 Jun 2024, 15:01
by mika
robtovey wrote: 04 Jun 2024, 11:00 So Material Pre Strain/Stress are only valid options when Eigen Analysis = True
Actually given Pre Stress/Strain affects the assembly in the second solution step also when

Geometric Stiffness = True

regardless of the value of the keyword Eigen Analysis. Then a prestressed stiffness matrix is created in the second step, but in the case of Eigen Analysis = False one must define some additional source term for the equations in order to get a nonzero solution (in your sif the given Strain Load serves for that purpose, so the second solution indeed changes when Pre Stress/Strain is altered). The eigen analysis doesn't need such source to obtain the eigenmodes corresponding to the prestressed stiffness matrix. It might be that the keyword Geometric Stiffness was originally added by having eigen analysis in mind.

If Geometric Stiffness = True, Pre Stress/Strain affects the assembly in the second solution step only, since the code has

Code: Select all

CALL StressCompose(..., GeometricStiffness .AND. iter>1,...)
Thus this option necessitates a two-step solution procedure.

Re: Stress Solver - Pre-strain is not applied?

Posted: 05 Jun 2024, 14:55
by robtovey
Hm, you're right, "Pre Strain" does have an effect on the simulations, but I can only see it when stress load is active too.

I've been playing with isotropic deformations of the unit disk because it's easy to track numerically and analytically (displacement of the perimeter is the same magnitude as the strain). I've attached the minorly updated script, you can just set the isotropic Pre Strain and Strain Load values in the top two lines.

The analytical equations are straight-forward. If I start with a disk of radius 1 and uniform isotropic strain, then I can rearrange

Code: Select all

 total_strain = (1 - analytical_radius) / analytical_radius 
to get

Code: Select all

 analytical_radius - 1 = - total_strain / (1 + total_strain) 
This passes my sanity checks (a disk of analytical_radius=1/2m stretched to a radius of 1m will have total_strain=1, similarly analytical_radius=2m compressed to 1m has total_strain = -0.5)

On the other hand, the numerical output I get in the VTU files is quite clear:

Code: Select all

simulated_radius - 1 = + strain_load / (1 + strain_load + pre_strain)
This equation is exact for the 12 pairs of (strain_load, pre_strain) that I've tried between -0.5 and 6.

Comparing the analytical and simulated radii/displacements, there are two things which don't look right to me:
  1. With pre_strain=0, the sign of the displacement and strain is wrong. The combination of strain_load = -0.5, pre_strain=0 predicts the disk shrinking to a point! The predicted strain is consistent with this (xx = yy = -1, xy = 0). If the sign had been flipped, then it would indicate exactly doubling the radius, which would be correct
  2. With large pre_strain, I would expect to get large displacement. Instead it seems to be missing from the top of the fraction (which is why it doesn't look like it's having an effect) which means it just dampens the strain load effect.
In both cases it looks like elmer is using "displacement = force / stiffness" where the stiffness matrix is being computed correctly but the force ignores Pre Strain and flips Strain Load.
Does this look like a bug to you?

Re: Stress Solver - Pre-strain is not applied?

Posted: 06 Jun 2024, 16:03
by mika
I also think that some signs in the code have not been consistent with the energy principle. As an attempt to correct this a fresh commit changes these signs as

https://github.com/ElmerCSC/elmerfem/co ... 2a46b846be

Re: Stress Solver - Pre-strain is not applied?

Posted: 07 Jun 2024, 10:47
by mika
If I run your case after these changes, with

Code: Select all

Pre Strain(3) = Real 5e-1 5e-1 0.0
...
Geometric Stiffness = True
...
Strain Load(3) = Real 5e-1 5e-1 0.0
I see a sensible result that the solution doesn't change in the two-step analysis

Code: Select all

ComputeChange: NS (ITER=1) (NRM,RELC): ( 0.25774787      2.0000000     ) :: mystress
ComputeChange: NS (ITER=2) (NRM,RELC): ( 0.25774787     0.17660338E-13 ) :: mystress
The first solution step now gives the displacement field u such that

stress = C E(u) - C E_0 = 0

where the initial strain E_0 corresponds to what has been given with the keyword Strain Load. In the second step the solution is done by adding the geometric stiffness matrix to the basic stiffness matrix. Since we have given

Code: Select all

Pre Strain(3) = Real 5e-1 5e-1 0.0
the linearization corresponding to the geometric stiffness matrix is done at the stress free-state, so the added geometric stiffness matrix is zero. That is, in the code the updates

Code: Select all

StressTensor(1,1) = StressTensor(1,1) + PreStress(1)
StressTensor(2,2) = StressTensor(2,2) + PreStress(2)
StressTensor(1,2) = StressTensor(1,2) + PreStress(3)
StressTensor(2,1) = StressTensor(2,1) + PreStress(3)
make the prestresses used in the computation of the geometric stiffness matrix zero. Therefore we basically solve the same system twice as also seen from the solver output.

The stress computation of the solver might need some rethinking so that the given initial strains and stresses would be taken into account.

Re: Stress Solver - Pre-strain is not applied?

Posted: 07 Jun 2024, 19:06
by robtovey
Thank you again for all your work on this.

I now get the updated equation

Code: Select all

simulated_radius - 1 = + strain_load / (1 + strain_load - pre_strain)
Honestly I don't really know what we're aiming for. What is your test to know if it's working correctly?
For example, without any documentation I don't know whether pre strain now has expected sign or not!

I still feel like a large pre strain with small strain load should correspond to a large translation, which is not currently the case, but I don't want this to drag on further than you do. If you think it's working fine, then feel free to say. Otherwise, I'm happy to continue to help debugging in any way I can.
Next week I will be on holiday though, so I'm sorry in advance for being slow to reply.