
[Sponsors] 
Support thread for "Solid Mechanics Solvers added to OpenFOAM Extend" 

LinkBack  Thread Tools  Display Modes 
April 28, 2014, 12:59 

#201  
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Quote:
I was aware of this issue with certain meshes, but I finally got around to looking at it and hopefully it is now solved with the attached detachFaceCracker.C. Try is out and let me know. Best regards, Philip 

April 30, 2014, 11:16 
Stress & Displacement Errors

#202 
New Member
Join Date: Mar 2014
Posts: 23
Rep Power: 4 
Hi,
A thin beam has been modeled by plane strain problem. When the number of cells along the thickness increases there is some error increasing in the simulation but by increasing number of cells along the length there is decreasing in the error, as shown in attached figures. Why does increasing cell numbers along the thickness cause an increasing in the error, please? Thank you. 

April 30, 2014, 11:44 

#203  
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Quote:
Using the finite volume method in its current linear form, large aspect ratio cells act stiffer in bending. So if you have long skinny cells along the length of the beam then it will act more stiff in bending than it should be. This behaviour can also occur in finite elements. I think someone has noted this earlier in this thread. Philip 

May 1, 2014, 22:08 

#204 
Member
Eric Bryant
Join Date: Sep 2013
Location: Texas
Posts: 44
Rep Power: 4 

May 9, 2014, 11:52 

#205  
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Quote:
You could try this case with a total displacement solver such as elasticNonLinTLSolidFoam instead of the incremental displacement solver and see how it compares. There are only a few changes to the case required: change all BCs with "nonLinear updatedLagrangian" to "nonLinear totalLagrangian"; and also change all "DU" to "U" in fvSchemes/fvSolution. Best regards, Philip 

May 10, 2014, 10:02 

#206 
New Member
Giampaolo Cetraro
Join Date: Oct 2012
Posts: 13
Rep Power: 5 
Hi Philip,
unfortunately there are no improvements using the total lagrangian solver. For the same timestep TL and UL give different behaviour. All the simulations use backward scheme and DU tolerance = 1e09. https://drive.google.com/file/d/0B1d...it?usp=sharing https://drive.google.com/file/d/0B1d...it?usp=sharing At least all the total lagrange simulations converge (with different modes) towards the steady state solution ...but now where the dissipation/dumping comes from? https://drive.google.com/file/d/0B1d...it?usp=sharing Best regards Giamp 

May 12, 2014, 04:16 

#207  
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Quote:
for the TL simulations using the backward scheme, did you modify the temporal term in the solver to be rho*d2dt2(U) instead of d2dt2(rho, U)? Philip 

May 12, 2014, 07:30 

#208 
New Member
Giampaolo Cetraro
Join Date: Oct 2012
Posts: 13
Rep Power: 5 
Hi Philip,
also with the temporal term rho*fvm::d2dt2(U) the solution is damped https://drive.google.com/file/d/0B1d...it?usp=sharing decreasing the time step, even if dumped, at least the solution seems to not change anymore https://drive.google.com/file/d/0B1d...it?usp=sharing but unfortunately at some some it crashes these are the last lines from the log of DeltaT=3.125e4 Code:
Time 1.87156, Corrector 9992, Solving for U using DICPCG, residual = 3.70566e05, relative residual = 2.45566, inner iterations 9 Time 1.87156, Corrector 9993, Solving for U using DICPCG, residual = 3.51757e05, relative residual = 0.244511, inner iterations 4 Time 1.87156, Corrector 9994, Solving for U using DICPCG, residual = 1.73622e05, relative residual = 1.30172, inner iterations 18 Time 1.87156, Corrector 9995, Solving for U using DICPCG, residual = 2.47067e05, relative residual = 0.183089, inner iterations 1 Time 1.87156, Corrector 9996, Solving for U using DICPCG, residual = 4.02368e05, relative residual = 0.197761, inner iterations 1 Time 1.87156, Corrector 9997, Solving for U using DICPCG, residual = 3.76838e05, relative residual = 0.389101, inner iterations 1 Time 1.87156, Corrector 9998, Solving for U using DICPCG, residual = 2.09545e05, relative residual = 0.672903, inner iterations 2 Time 1.87156, Corrector 9999, Solving for U using DICPCG, residual = 1.19139e05, relative residual = 0.731561, inner iterations 6 Time 1.87156, Solving for U, Initial residual = 0.0049075, Final residual = 1.19139e05, Relative residual = 0.731561, No outer iterations 10000 ExecutionTime = 2376.78 s ClockTime = 2398 s ExecutionTime = 2376.78 s Time: 1.87187 Time 1.87187, Corrector 0, Solving for U using DICPCG, residual = 0.00534617, relative residual = 1, inner iterations 1 Time 1.87187, Corrector 1, Solving for U using DICPCG, residual = 1.96508e05, relative residual = 1.01578, inner iterations 14 Time 1.87187, Corrector 2, Solving for U using DICPCG, residual = 2.91777e05, relative residual = 5.95844, inner iterations 5 Time 1.87187, Corrector 3, Solving for U using DICPCG, residual = 2.74333e05, relative residual = 1.4765, inner iterations 2 Time 1.87187, Corrector 4, Solving for U using DICPCG, residual = 2.01264e05, relative residual = 0.445946, inner iterations 1 Time 1.87187, Corrector 5, Solving for U using DICPCG, residual = 1.45802e05, relative residual = 4.99339, inner iterations 4 Time 1.87187, Corrector 6, Solving for U using DICPCG, residual = 2.34268e05, relative residual = 0.336645, inner iterations 1 Time 1.87187, Corrector 7, Solving for U using DICPCG, residual = 2.84046e05, relative residual = 0.698792, inner iterations 2 Time 1.87187, Corrector 8, Solving for U using DICPCG, residual = 3.61341e05, relative residual = 2.96854, inner iterations 5 Time 1.87187, Corrector 9, Solving for U using DICPCG, residual = 3.78643e05, relative residual = 2.09039, inner iterations 2 Time 1.87187, Corrector 10, Solving for U using DICPCG, residual = 3.68786e05, relative residual = 4.16278, inner iterations 7 Time 1.87187, Corrector 11, Solving for U using DICPCG, residual = 3.29223e05, relative residual = 1.70304, inner iterations 4 Time 1.87187, Corrector 12, Solving for U using DICPCG, residual = 3.08959e05, relative residual = 9.57103, inner iterations 8 Time 1.87187, Corrector 13, Solving for U using DICPCG, residual = 4.94034e05, relative residual = 0.192556, inner iterations 1 Time 1.87187, Corrector 14, Solving for U using DICPCG, residual = 9.53101e05, relative residual = 0.190417, inner iterations 1 Time 1.87187, Corrector 15, Solving for U using DICPCG, residual = 0.000259233, relative residual = 0.18017, inner iterations 1 Time 1.87187, Corrector 16, Solving for U using DICPCG, residual = 0.00113276, relative residual = 1.74945, inner iterations 1 Time 1.87187, Corrector 17, Solving for U using DICPCG, residual = 0.0123064, relative residual = 1.9589, inner iterations 1 Time 1.87187, Corrector 18, Solving for U using DICPCG, residual = 0.377465, relative residual = 1.48288, inner iterations 1 Time 1.87187, Corrector 19, Solving for U using DICPCG, residual = 0.995577, relative residual = 1.00002, inner iterations 1 Time 1.87187, Corrector 20, Solving for U using DICPCG, residual = 1, relative residual = 1, inner iterations 1 > FOAM Warning : From function eigenValues(const tensor&) in file primitives/Tensor/tensor/tensor.C at line 170 complex eigenvalues detected for tensor: (5.24059e+18 1.52383e+20 6.46432e26 378752 1 4.87694e11 653.646 18991.2 1) > FOAM Warning : From function eigenValues(const tensor&) in file primitives/Tensor/tensor/tensor.C at line 170 complex eigenvalues detected for tensor: (5.24059e+18 1.52383e+20 6.46432e26 378752 1 4.87694e11 653.646 18991.2 1) https://drive.google.com/file/d/0B1d...it?usp=sharing Best regards, Giamp 

May 12, 2014, 07:36 

#209 
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Hi Giamp,
Hmnnn for transient analysis the beam should never reach steadystate, can you check if a small bugfix in TL solver makes a difference: in writeFields.H, change this: Code:
rho = rho/J; Code:
//rho = rho/J; Philip 

May 12, 2014, 11:50 

#210 
New Member
Giampaolo Cetraro
Join Date: Oct 2012
Posts: 13
Rep Power: 5 
good! with the bugfix, the beam keeps oscillating
https://drive.google.com/file/d/0B1d...it?usp=sharing but still, decreasing the timestep, at some point the solutions are different. ElasticSolidFoam keeps oscillating but crashes when decreasing the timestep https://drive.google.com/file/d/0B1d...it?usp=sharing 

May 13, 2014, 22:11 

#211 
Senior Member
Daniel
Join Date: Mar 2013
Posts: 211
Rep Power: 13 
Hi Dr. Cardiff
Have you ever tried the FSI solver using velocityLaplacian mesh motion solver in parallel? I am facing an error: Code:
Constructing processor meshes Processor 0 Number of cells = 10470 Number of faces shared with processor 1 = 171 Number of processor patches = 1 Number of processor faces = 171 Number of boundary faces = 21405 Processor 1 Number of cells = 10470 Number of faces shared with processor 0 = 171 Number of processor patches = 1 Number of processor faces = 171 Number of boundary faces = 21395 Number of processor faces = 171 Max number of processor patches = 1 Max number of faces between processors = 171 Processor 0: field transfer > FOAM FATAL ERROR: size of field = 21744 is not the same as the size of mesh = 21578 From function DimensionedField<Type, GeoMesh>::DimensionedField(const IOobject& io,const Mesh& mesh, const dimensionSet& dims, const Field<Type>& field) in file /opt/OpenFOAM/OpenFOAM3.0ext/src/OpenFOAM/lnInclude/DimensionedField.C at line 71. FOAM aborting Aborted (core dumped) Best wishes 

May 16, 2014, 14:00 

#212 
New Member
Giampaolo Cetraro
Join Date: Oct 2012
Posts: 13
Rep Power: 5 
Hi Philip,
I'm trying to modify the solver elasticNonULSolidFoam using an incompressible MooneyRivlin hyperelastic model. From the strain energy function W = c1*(tr(C)  3) + c2*(0.5*((tr(C))^2  tr(C^2)) 3) I derived the second PiolaKirchhoff stress tensor in its incremental form DSigma = 2*c1*I + 2*c2*tr(C)*I  2*c2*C  p*inv(C) where C is the right CauchyGreen deformation tensor (C=F.T() & F), and c1+c2=mu/2 Following http://powerlab.fsb.hr/ped/kturbo/Op...tressPaper.pdf and http://powerlab.fsb.hr/ped/kturbo/op...18062007.pdf I came up with the following Updated lagrangian momentum equation (neglecting the pressure term) rho*fvm::d2dt2(DU) ==  fvm::laplacian(2*c2, DU, "laplacian(DDU,DU)") + fvc::div ( 4*c2*tr(gradDU)*I + 2*c2*(gradDU.T()), "div(sigma)" ) + fvc::div ( 2*c2*(I*tr(gradDU.T() & gradDU))  2*c2*(gradDU.T() & gradDU) + ((sigma + DSigma) & gradDU), "div(sigma)" ) Basically in elasticNonLinULSolidFoam I changed only the momentum equation and DSigma (neglecting the pressure term for now). I calculated C as C = (2*DEpsilon + I) and used c1=0.4375*mu and c2=0.0625*mu. I tried to run a steady state simulation with this solver with a case which works fine with elasticNonLinULSolidFoam, but unfortunately the simulation crashes Code:
Time 0.0001, Corrector 61, Solving for DU using GAMG, res = 0.786476, rel res = 0.117896, inner iters 1 Time 0.0001, Corrector 62, Solving for DU using GAMG, res = 0.783442, rel res = 0.116527, inner iters 1 Time 0.0001, Corrector 63, Solving for DU using GAMG, res = 0.780196, rel res = 0.115416, inner iters 1 Time 0.0001, Corrector 64, Solving for DU using GAMG, res = 0.776698, rel res = 0.116016, inner iters 1 Time 0.0001, Corrector 65, Solving for DU using GAMG, res = 0.77302, rel res = 0.156058, inner iters 1 Time 0.0001, Corrector 66, Solving for DU using GAMG, res = 0.799312, rel res = 0.628151, inner iters 1 Time 0.0001, Corrector 67, Solving for DU using GAMG, res = 0.998106, rel res = 0.998607, inner iters 2 Time 0.0001, Corrector 68, Solving for DU using GAMG, res = 1, rel res = 0.999996, inner iters 2 Time 0.0001, Corrector 69, Solving for DU using GAMG, res = 1, rel res = 1, inner iters 2 Time 0.0001, Corrector 70, Solving for DU using GAMG, res = 1, rel res = 1, inner iters 2 Floating point exception (core dumped) Any suggestion to make it work? About the pressure term, since that should be an arbitrary pressure, does assuming p=0 make sense? ...or the problem might come from this assumption? Thank you very much! 

May 20, 2014, 13:12 

#213 
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Hi Giampaolo,
See this paper by Bijelonja et al, they develop the mathematical model for an incompressible MooneyRivlin material using a large strain total Lagrangian approach. They use the SIMPLE algorithm to couple the pressure term, just like in fluids. So essentially merge elasticNonLinTLSolidFoam and simpleFoam and change calculation of stress in terms of MooneyRivlin. Also traction boundary gradient method must be changed. Let me know how it goes, I previously looked at creating a solid solver for incompressible materials using SIMPLE/PISO coupling and I had problems with traction boundary conditions, but it is possible as the authors above have achieved it using cellcentred Finite Volume method. Good luck, By the way, as regards the timestep defence of the oscillating membrane problem, I am not sure what the issue is, I suppose a more accurate time discretisation (e.g. Newmark, CrankNicholson, etc.) would help but I am not sure. Philip 

May 23, 2014, 13:23 
MesquiteMotionSolver with icoFsiElasticNonLinULSolidFoam

#214 
New Member
Join Date: Feb 2014
Posts: 10
Rep Power: 4 
@Philip
I wanted to remesh the deformed part of the mesh while using icoFsiElasticNonLinULSolidFoam but when i used dynamicTopoFvMesh along with mesquiteMotionSolver I got this error msg!!! > FOAM FATAL ERROR: Problem with mesh motion solver selection From function icoFsiElasticNonLinULSolidFoam in file moveFluidMesh.H at line 140. FOAM aborting Aborted (core dumped) and when I forced the moveFluidMesh.H to use the given motion solver, it reported an error msg.: > FOAM FATAL ERROR: request for tetPointVectorField motionU from objectRegistry region0 failed available objects of type tetPointVectorField are 0 ( ) From function objectRegistry::lookupObject<Type>(const word&) const in file /home/a/foam/foamextend3.0/src/foam/lnInclude/objectRegistryTemplates.C at line 139. FOAM aborting can I know where the problem is? and what should I do to solve this issue? regards, MA Last edited by sfmoabdu; May 25, 2014 at 03:45. 

May 28, 2014, 08:05 
Dual BC

#215 
New Member
Join Date: Mar 2014
Posts: 23
Rep Power: 4 
Hi,
How can a traction and displacement BCs be applied simultaneously in a face, please? For example, having tractionDisplacement and slip types BC in the same time in a face. Thanks. 

May 28, 2014, 08:54 

#216  
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Quote:
Philip 

May 28, 2014, 10:15 
Dual BC

#217 
New Member
Join Date: Mar 2014
Posts: 23
Rep Power: 4 
Actually there are shear and displacement in a single face, as shown in the left BC in attachment.


May 28, 2014, 10:20 

#218  
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 
Quote:
there are 'frictionless' rollers on the left which implies there is zero shear? And there is a fixed point at the centre of the left boundary. If this is the case, then you would split the left boundary so as there is a thin fixedDisplacement boundary for the fixed part and fixedDisplacementZeroShear for the other parts. Philip Last edited by bigphil; May 28, 2014 at 10:21. Reason: typo 

May 28, 2014, 10:31 
Dual BC

#219 
New Member
Join Date: Mar 2014
Posts: 23
Rep Power: 4 
Please look at the paper in attachment. I would like to know how beams D and E in Fig. 2 on page 2 can be modeled by stressedFoam. Applying the left dual BCs is my question.
Thanks. 

May 28, 2014, 10:34 

#220 
Senior Member
Philip Cardiff
Join Date: Mar 2009
Location: Dublin,Ireland
Posts: 617
Rep Power: 21 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
ESIOpenCFD Releases OpenFOAM v2.2.0  opencfd  OpenFOAM Announcements from ESIOpenCFD  13  March 30, 2013 17:52 
[Virtualization] OpenFOAM oriented tutorial on using VMware Player  support thread  wyldckat  OpenFOAM Installation  2  July 11, 2012 16:01 
GPU Linear Solvers for OpenFOAM  gocarts  OpenFOAM Announcements from Other Sources  35  March 1, 2012 21:41 
Crosscompiling OpenFOAM 1.7.0 on Linux for Windows 32 and 64bits with Mingww64  wyldckat  OpenFOAM Announcements from Other Sources  3  September 8, 2010 06:25 
OpenFOAM Debian packaging current status problems and TODOs  oseen  OpenFOAM Installation  9  August 26, 2007 13:50 