Quote:
Thank you for providing me with the information and the links. By mixed I mean that it is a velocity/strains formulation where the primary variables are linear momentum and the deformation gradient tensor. Hopefully, I will be able to share more details regarding our method in this thread very soon. Cheers! |
crackStressFoam
Hello there,
(forgot to say its the crackStressFoam in extend-bazaar FluidStructureInteraction package) I'm interested in the crackStressFoam solver (especially for FSI) but I cannot find out anything about it anywhere .. apart from the .C file - Finite volume stress analysis solver for arbitrary crack propagation Does anybody have any links or tutorials? As it resides in the FluidStructureInteraction directory I'm assuming I can use it for FSI? Any advice appreciated. Br Jason |
thermalStressFoam hotCylinder validation
2 Attachment(s)
Dear Phil,
It would seem that people succeed in validating the hotCylinder case with the FSI thermalStressFoam solver. Your plot here among others. I've not yet managed to validate thermalStressFoam with the hotCylinder case. Bug fix done. I have uploaded my case for scrutiny. Any feedback would be greatly appreciated. Kind regards, Will. |
Default thermalStressFoam hotCylinder validation II
1 Attachment(s)
The problem I had was one of interpretation.
In the case I uploaded last week the cylinder is not free to deform longitudinally (2-D plane strain), whereas in the case of Timoshenko the cylinder is free to deform longitudinally. Cheers, Will. |
fsiFoam: Restart-Issue
Hi,
has anyone successfully solved the restart-issue using fsiFoam ? I tried to restart beamInCrossFlow using foam-extend-3.1 and fsiFoam from extend-bazaar. It does restart at several different time steps, but looking at the new time steps after restart, one can see a "jump" in the displacement of the fluid-interface and also a "gap" between the fluid-interface and the solid-interface. Also the solution at the end of the run (static situation) is different from the original behaviour (without restart). Should that work in theory ? Has anyone similar problems ? Thanx for help in advance, Martin Quote:
|
1 Attachment(s)
Dear all
has anyone faced the problem of interface separation? , in which the solid interface moves away from fluid interface (as in the attached snapshot) ... it happens very slowly without any crash of solution and at some time step the gap becomes very big and the solver gives "Floating point exception (core dumped) " or "complex eingenvalues in the matrix ... " !! has anyone of you the solution of that? it happens in both fsiFoam and icoFsiElasticNonLinULSolidfoam. by the way i have used the HronTurek tutorial but i just changed the material and thickness of the plate. Furthermore, the fluid domain is perfectly solved as only CFD. |
New article on large strain computational solid dynamics
Dear All,
An article on large strain computational solid dynamics has recently been published in International Journal for Numerical Methods in Engineering. Please find below a link to the manuscript on researchgate. https://www.researchgate.net/publica...train_dynamics Thanks |
Quote:
The procedures look very interesting. Do you intend to share the code? Was it implemented in OpenFOAM? - from the article - yes All the best, Philip |
Quote:
Hi Philip, Thanks for the appreciation. Yes, you are correct. The methodology has been implemented from scratch in OpenFOAM. Currently, we are working on the code to extend it to contact mechanics applications. Once that is achieved we plan to share the code. |
Restart FsiFoam
Hi Philip,
did you have time to have a look at the restart-issue ? Would be really nice if it would work.. Martin Quote:
|
Quote:
|
Quote:
Regards, Maimouna |
Quote:
Hi Maimouna, here's my case, the data of the geometry comes from a benchmark case, but the solid part deforms little. Thanks, |
Dear Phillip,
I am using the elasticSolidFoam (foam extended 3.1) solver to apply traction to a "simple" system. With simple, I mean geometrically because it is just a plate that's fixed at the bottom and pulled at the top. This plate though has residual stresses within, so a "body force". I added to the momentum equation this body force that is in the material as an extra term (rho*div(residual stress)). I read in one of the forums for a different case, where the gravity was supposed to be introduced, that you said this force is added only in the first time step, nonetheless I just included directly in the equation (mistake??). It looks like: fvVectorMatrix UEqn ( rho*fvm::d2dt2(U) == fvm::laplacian(2*muf + lambdaf, U, "laplacian(DU,U)") + divSigmaExp + rho*divStress (Body force) ); What I find so far is that the Young modulus of my system changes when I change the fixedDisplacement (strain on top) when this body force in present. When it is not present it doesn't change. So there is a dependence on the strain (rate). My first question is if adding the body force like I did in the momentum equation makes sense?? :confused: Also I am puzzled with this strain dependence of the elastic modulus because there is no time scale in the momentum equations. :confused: It would be very nice hearing from you, Kind regards, Heliana |
Quote:
Is your residual stress a uniform field or do you have some spatially varying data? To add a residual stress, you would add the divergence of this stress to the momentum equation: Code:
fvVectorMatrix UEqn Update the definition of stress in the "elasticSolidFoam/calculateEpsilonSigma.H": Code:
epsilon = symm(gradU); Code:
// Add residual stress term to boundary gradient As regards applying the term in just the first time-step, this would only be the case for an incremental solver where the momentum equation is cast in incremental form; for elasticSolidFoam it solves the conservation of total linear momentum so the term should be added every time-step. For the time dependence, there is a temporal term in the momentum equation 'd2dt2(U)'; this term may be significant for highly transient problems or over a very short time. The implemented constitutive law (Hooke's law) is rate independent. Also, you should check that the solver is converging in each time-step. Philip |
Hi Philip,
Thank you very much for answering. My field is not uniform across the geometry but varies spatially. I am implementing the changes you suggested, I think I was confused with the other forum about the incremental solver, so now I understand that the "myResidualStressField" should be added every time-step. As soon as it works I let you know! I struggle in the moment though when compiling the $FOAM_SRC/solidModels/constitutiveModel/tractionBoundaryGradient/tractionBoundaryGradient.C I get this error: constitutiveModel/tractionBoundaryGradient/tractionBoundaryGradient.C: In static member function ‘static Foam::tmp<Foam::Field<Foam::Vector<double> > > Foam::tractionBoundaryGradient::snGrad(const vectorField&, const scalarField&, const Foam::word&, const Foam::word&, const Foam::fvPatch&, bool, const Foam::nonLinearGeometry::nonLinearType&, bool)’: constitutiveModel/tractionBoundaryGradient/tractionBoundaryGradient.C:541:17: error: invalid initialization of reference of type ‘const volSymmTensorField& {aka const Foam::GeometricField<Foam::SymmTensor<double>, Foam::fvPatchField, Foam::volMesh>&}’ from expression of type ‘bool’ ( Anyways, I will find a way to work it out I let you know how he solver converges. Cheers, Heliana |
Hi Heliana,
For the compilation error, this was my fault, there was a typo in the code I posted: This line: Code:
patch.boundaryMesh().mesh().foundObject<volSymmTensorField> Code:
patch.boundaryMesh().mesh().lookupObject<volSymmTensorField> Philip |
Hi Philip,
I fixed it like you said and I compiled and runs fine. This is just a bit of brainstorming: I was thinking though that at t=0 the only stresses present would be my residual stresses which must then create an initial deformation -- sort of saying it comes inherently within the material -- so I figure I have to initialize or calculate somehow that initial deformation where the calculation must start from at the beginning. I think what I mean is that the initial deformation, when there are residual stresses present cannot be zero but a value that should be calculated with the elasticSolidFoam solver assuming that sigma is zero by myResStressField is not zero, hence I should initialize or calculate this deformation a the beginning right and then start from there right? Sorry if I put it too complicated :) Cheers, Heliana |
Sorry for the messed message,
After thinking about it I think I can rephrase it better. Since I have residual stresses in my geometry, I could calculate an initial deformation assuming linear elasticity - Hooke's Law - and this initial deformation would be the initial internal field I would use in U to start my calculation which also includes the residual stresses in the momentum equation. So I would calculate an initial deformation called DU, hence: volTensorField gradDU = fvc::grad(DU); // Hooke's law myResStress = 2*mu*(symm(gradDU)) + lambda*(I*tr(symm(gradDU))) I need somehow to calculate DU from the previous equation and then that would be my initial inherent deformation across my geometry (because of those residual stresses), so then it is like my system is initialized with this field and then I use the solver.... Does it make sense at all? or am I just running in circles? Heliana |
Quote:
Yes, I understand; if the residual stress field is not in equilibrium then it would cause some deformation. With the suggested modifications, the elasticSolidFoam solver will calculate these deformations in the first time-step. For example, to see this, you could set all your boundary conditions to traction-free (i.e. solidTraction with all zeros) and then run the solver for one time-step; you may need to set one patch to fixedDisplacement to stop the model floating away. If the residual stress field was not in equilibrium then you will see that the model will deform to reach an equilibrium between the residual stresses and generated internal elastic stresses. In short, I don't think you need to make any additional changes, the solver already solves for and accounts for this residual deformation. Philip |
Hi Philip,
Yes I also thought something like that but I was getting a Floating point exception error. I am going to work on it and keep you updated :) Cheers, Heliana |
Hi all,
FYI, a couple of recent papers on solid mechanics in OpenFOAM:
I am in the process of tidying things up and intend to publicly share the main procedures. PDFs of the accepted articles (pre publishing) can be found at my ResearchGate: https://www.researchgate.net/profile/Philip_Cardiff Philip |
geomechanics|approach to compute body forces/weight as initial cond.
Hello out there,
this seems to be a space where the solid mechanics fraction seems to interchange. However I am not sure if this is the right place. I'd like to apologize if this isn't and I'd appreciate if somebody told me where to post. I have been exploring Tian Tangs implementation for geomechanical problems, which has proved easier than expected due to the excellent documentation. The analysis concerns an (almost) water saturated porous matrix. Coming from geohydraulics I have to adapt myself to the incremental formulation used in the solid mechanic solvers, so my question might sound a bit silly: I wish to assign realistic initial conditions to pore pressure p and stress sigma which stem from gravity. However in the incremental formulation of the solver body forces are neglected. I have no clue as how to assign water and soil a self-weight in order to obtain a realistic initial condition, hydrostatic for the fluid and geostatic for the soil. This seems necessary for the analysis of plastic deformations, which response to the total and not only to the incremental stress changes. I am trying to reproduce the wave-seabed-gravity structure interaction tutorial in which a consolidation is run to provide initial conditions applying a "K0-procedure", which unfortunately is not explained further. Must I apply some kind of operation with the gravity vector in order to assign body forces to matrix and fluid? How should such an approach be implemented? An other option could be to take the z coordinates from constant/polyMesh and try to generate hydrostatic and geostatic fields in octave/excel but his approach seems not too appealing. I could imagine, that this question has been raised by other so any hints appreciated on the "K0 approach" or any other way to compute body forces/weight in order to start the analysis from realistic initial conditions. Tom |
Quote:
It's possible to include an initial stress and initial pore-pressure by writing the conservation equation in total (not incremental) form, but you can still solve for the increment of displacement for convenience i.e. div(totalStress) == 0 div(initialStress + effectiveStress + initialPorePressure + porePressure) == 0 div(initialStress + oldTime_effectiveStress + increment_effectiveStress + initialPorePressure + porePressure) == 0 so for example, in OpenFOAM this could be: Code:
fvVectorMatrix DUEqn Philip |
Quote:
sorry for late reply, if you still looking for help email me please may78may@hotmail.com. Best, Maimouna |
How is traction from fluid applied on solid??
Hello all,
I am using the icoFsiElasticNonLinULSolidFoam solver. I cannot understand how the traction force (tForce) calculated in setInterfaceForce.H is applied in the solveSolid.H. The solveSolidEuler.H has the governing equation implemented, but i dont see the force term in the governing equation. Can anyone please explain the steps in which solid is solved? How is that traction force applied? and where? Please guide me on this. Following code is solveSolidEuler.H { # include "readSolidMechanicsControls.H" int iCorr = 0; lduMatrix::solverPerformance solverPerf; scalar initialResidual = 0; lduMatrix::debug = 0; # include "EulerCoeffs.H" do { DU.storePrevIter(); fvVectorMatrix DUEqn ( Cn*rho*fvm::ddt(DU) - Co*rho*DV.oldTime() == fvm::laplacian(2*mu + lambda, DU, "laplacian(DDU,DU)") - fvc::laplacian(mu + lambda, DU, "laplacian(DDU,DU)") + fvc::div ( mu*gradDU.T() + lambda*(I*tr(gradDU)) + mu*(gradDU&gradDU.T()) + 0.5*lambda*(I*tr(gradDU & gradDU.T())) + (sigma & DF.T()) + (DSigma & DF.T()), "div(sigma)" ) ); solverPerf = DUEqn.solve(); DU.relax(); if(iCorr == 0) { initialResidual = solverPerf.initialResidual(); } gradDU = fvc::grad(DU); DF = gradDU.T(); # include "calculateDSigma.H" } while ( solverPerf.initialResidual() > convergenceTolerance && ++iCorr < nCorr ); Info << "Solving for " << DU.name() << ", Initial residual = " << initialResidual << ", Final residual = " << solverPerf.initialResidual() << ", No outer iterations " << iCorr << endl; DV = fvc::ddt(DU); lduMatrix::debug = 1; } |
Thank you, Philip.
That approach is exactly what I was looking for. I will try to implement that (which might take me a while due to my limited oF experience) and will report the (hopefully successful) result. |
Quote:
For example: If you installed OpenFOAM- 2.3.0 and foam- extend 3.1 in your machine, you have to start the terminal with of230 --> to start using OpenFOAM- 2.3.0 OR fe31 --> to start foam-extend 3.1. Kind regards Maimouna |
Hello everyone,
Did anybody have any success with the restarting issue for FSI cases running in parallel, discussed before on this thread? I am currently using foam-extend-4.0 with the bazaar extension to simulate FSI problems, but I am still facing the same errors requiring the V0 file for the solid part (it happened too with the last 3.2 version). I tried the workaround proposed by Needled, but still did not work. Thanks in advance :) Iago |
Quote:
I tried to run some cases using fsiFoam solver in parallel, but that wouldn't help me much. I mean in speed up. Also, simulation in parallel gives different results when compared with serial simulation, because of the number of cores. Briefly, I think, foam-extend 4.0 doesn't implement to work in parallel, it still needs some modification. Kind regards Maimouna |
Quote:
But were you able to restart a case decomposed in fsiFoam? Always when I try it, I get the error: Code:
[7] Regarding the speedup, with all the cases I worked with, it was better to run in parallel. But thanks for pointing out that you were having different results with serial and parallel, I need to compare my results to check it! Regards, Iago |
Hello,
There was a question about dynamicTopoFvMesh, but it was left without an answer. I'm working on a case similar to "beamInCrossFlow", but with a more complicated geometry, so I need to add to it the posibility of making topological changes to the fluid mesh. I had some troubles with the parameters, "cacheAgglomeration" in fvSolution, etc., but finally I've stucked on a problem with GGI. Each time after remeshing I'm getting an error: Quote:
I'm using "velocityLaplacian" motion solver, since "mesquite" seems not to be availble with fsiFoam ("Problem with fluid mesh motion solver selection"). Is there any possibility to run fsiFoam with a changing topology? I would be very grateful for any clues. Kind regards, Bartosz Żłobiński |
Quote:
Have you tried deleting the V0 file from the latest fluid folder in all processor directories? I've found this to work for me. I think you may need to delete the meshPhi files from the latest fluid folders, too. KR |
Hi all,
I've been trying to implement LES to the "beamInCrossFlow" tutorial for the FSI toolkit for foam-extend 4.0. My problem is that a heavily oscillating pressure field appears in the fluid field, eventually causing extreme velocities, causing the simulation to crash. My changes with respect to the original tutorial are as follows: - Rewritten blockMeshDict to a complete 3D domain for both the fluid and structure. - Increased the max speed in setInletVelocity to 10, in order to get turbulence. - Changed the properties of the solid to the following: rho = 8000, E = 4e6, nu = 0.3 - Implemented PISO solver with LES model oneEqEddy, with vanDriest smoothing. - Decreased FSI tolerance to 1e-5 and increased nOuterCorr to 200. - ddtScheme: CrankNicolson 0.8 - div(phi,U) Gauss filteredLinear2 1 0 - fvSolution changed to common settings for LES simulations. - The Courant number has been kept at 0.4-0.5. The ddtScheme was chosen by trying out all the ddtSchemes without FSI. The use of other ddtSchemes, including CrankNicolson 0.7 or 0.9, causes the simulations without FSI to blow up, too. Simulations with CrankNicolson 0.8 run smoothly without FSI however. In all cases, the problem has been that a field of high pressure at the lower edge of the inlet has merged with the field of high pressure in front of the beam. When these two fields of high pressure merge, they start to "pump" the fluid over the beam, resulting in pressure oscillations for the entire field. The pressure typically ends up oscillating from +200 Pa to -200 Pa for 80% of the internal field within 0.05 seconds. I therefore assumed that the inlet was too close to the beam, and increased the length from the inlet by 0.5m. The pressure oscillations did still occur. I have also tried the following adjustments: - Max speed 5 m/s. - Initial internal fields calculated by consistentIcoFluid with FSI, by LES simulation without FSI, and by increasing the velocity from zero with LES and FSI. - Reduced cell lengths to the half in all directions. - Adjusted time step to Co = 0.2 and Co = 0.8. In all these configurations, the pressure oscillations do still occur. The pressure oscillations do not occur when doing a LES simulation of the undeformed geometry without FSI, and neither if the stiffness (E) of the solid is set to be very high in the FSI simulation (causing it to be stationary). Has anyone else tried to implement LES to this case? I would be very grateful for any suggestions. Best regards, KR |
Quote:
The problem is with the globalFaceZones: after a topological change, these get messed up and it is not trivially to fix them; consequently, as you have seen the GGI will complain. In your case, is there a topo change at the patch/faceZone? If so, have you tried a case when the topo change is away from the patch? Just an idea. Similar issue discussed here: https://sourceforge.net/p/openfoam-e...ndrelease/280/ Philip |
Thank you for your explanation and suggestion, bigphil. I'll try to modify the dynamicMeshDict to force the topo changes to happen only far enough from the interface (e.g. two rows of cells), it should be sufficient for my simulation.
Hi, rege. I was trying to use LES for my case similar to the "beamInCrossFlow" and it seems to work fine. I have also encountered extreme velocities, but they haven't been caused by such pressure oscilations and apeared also for DNS. Kind regards, Bartosz Żłobiński |
Hello again,
Unfortunately, I still can't make my case working with dynamicTopoFvMesh. I've switched off edgeRefinement and only some swapping occurs on the inlet, outlet and on the internal cells. I'm still getting the same error from GGI, although there are no changes on the interface patch (it's named "granica"). I'm enclosing my dynamicMeshDict. Code:
/*--------------------------------*- C++ -*----------------------------------*\ Kind regards, Bartosz Żłobiński |
Quote:
Thank you very much for the answer! I tried deleting the all the old time-steps information files from the last time-step folder and didn't work. But deleting them and also the meshPhi file finally worked :) Regards, Iago |
problem with solverPerformance-class in porous media flow/deformation
Hi outthere,
due to an upgrade (opensuse 13.2 -> leap 42.2) I had to reinstall openFoam. Being working with oF3.0ext I decided to got for oF4.0ext. After installation everything seemed to work well and I could test a couple of tutorials (icoFoam, plateHole etc). However, when I tried to carry on with modifications for elastoPlastoBiotFoam (coupled flow and deformation in prorous media) some compiling errors appeared, not seen before: Code:
pDependElastPlastBiotFoam.C: In function ‘int main(int, char**)’: I have no idea, why this was never a problem in oF3.0ext and have no clue what is going out. I suppose i should modify the MAKE/options, which actually are set to: Code:
EXE_INC = \ Tom |
Quote:
In foam-extend-4.0, these lines: Code:
lduMatrix::solverPerformance solverPerfP; Code:
lduSolverPerformance solverPerfP; Code:
word correct=word(runTime.controlDic t().lookup("correctPlasticity")); Philip |
All times are GMT -4. The time now is 04:18. |