From a3b9ee8fac50db3decf0715003e3f05dcb819923 Mon Sep 17 00:00:00 2001 From: Hrvoje Jasak Date: Thu, 10 Apr 2014 17:45:26 +0100 Subject: [PATCH] Update and clean-up of solid mechanics --- applications/solvers/solidMechanics/Allwmake | 2 + .../solidMechanics/deprecatedSolvers/Allwmake | 2 - .../stressFemFoam/calculateStress.H | 53 ------ .../calculateDEpsilonDSigma.H | 4 +- .../calculateDivDSigmaExp.H | 101 +++++----- .../calculateDivDSigmaNonLinExp.H | 103 ++++++---- .../elasticPlasticNonLinTLSolidFoam.C | 178 +++++++++--------- .../readDivDSigmaExpMethod.H | 25 ++- .../writeFields.H | 159 ++++++++-------- .../stressFemFoam/Make/files | 0 .../stressFemFoam/Make/options | 0 .../stressFemFoam/calculateStress.H | 120 ++++++++++++ .../stressFemFoam/createFemFields.H | 0 .../stressFemFoam/createFemMesh.H | 0 .../stressFemFoam/createFields.H | 0 .../Traction/TractionPointPatchVectorField.C | 13 +- .../Traction/TractionPointPatchVectorField.H | 0 .../tractionTetPolyPatchVectorField.C | 0 .../tractionTetPolyPatchVectorField.H | 0 .../stressFemFoam/readControls.H | 0 .../stressFemFoam/readMechanicalProperties.H | 0 .../stressFemFoam/stressFemFoam.C | 4 +- .../stressFemFoam/Allclean | 0 .../stressFemFoam/Allrun | 0 .../stressFemFoam/plateHole/0/U | 0 .../plateHole/constant/mechanicalProperties | 0 .../plateHole/constant/polyMesh/blockMeshDict | 60 +++--- .../plateHole/constant/polyMesh/boundary | 2 +- .../plateHole/system/controlDict | 6 +- .../stressFemFoam/plateHole/system/fvSchemes | 0 .../stressFemFoam/plateHole/system/fvSolution | 0 .../plateHole/system/tetFemSolution | 10 +- 32 files changed, 483 insertions(+), 359 deletions(-) delete mode 100644 applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/calculateStress.H rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/Make/files (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/Make/options (100%) create mode 100644 applications/solvers/solidMechanics/stressFemFoam/calculateStress.H rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/createFemFields.H (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/createFemMesh.H (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/createFields.H (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.C (96%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.H (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/femStress/tractionTetPolyPatchVectorField.C (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/femStress/tractionTetPolyPatchVectorField.H (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/readControls.H (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/readMechanicalProperties.H (100%) rename applications/solvers/solidMechanics/{deprecatedSolvers => }/stressFemFoam/stressFemFoam.C (97%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/Allclean (100%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/Allrun (100%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/0/U (100%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/constant/mechanicalProperties (100%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/constant/polyMesh/blockMeshDict (72%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/constant/polyMesh/boundary (96%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/system/controlDict (94%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/system/fvSchemes (100%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/system/fvSolution (100%) rename tutorials/solidMechanics/{deprecatedTutorials => }/stressFemFoam/plateHole/system/tetFemSolution (86%) diff --git a/applications/solvers/solidMechanics/Allwmake b/applications/solvers/solidMechanics/Allwmake index 7c7dc0f32..58381dd5b 100755 --- a/applications/solvers/solidMechanics/Allwmake +++ b/applications/solvers/solidMechanics/Allwmake @@ -20,5 +20,7 @@ wmake elasticThermalSolidFoam wmake icoFsiElasticNonLinULSolidFoam wmake viscoElasticSolidFoam +wmake stressFemFoam + (cd utilities ; wmake all) (cd deprecatedSolvers ; ./Allwmake) \ No newline at end of file diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/Allwmake b/applications/solvers/solidMechanics/deprecatedSolvers/Allwmake index 1e2a5ecab..612815a04 100755 --- a/applications/solvers/solidMechanics/deprecatedSolvers/Allwmake +++ b/applications/solvers/solidMechanics/deprecatedSolvers/Allwmake @@ -11,8 +11,6 @@ wmake contactStressFoam wmake newStressedFoam wmake newContactStressFoam -wmake stressFemFoam - wmake icoFsiFoam wmake solidDisplacementFoam diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/calculateStress.H b/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/calculateStress.H deleted file mode 100644 index f81f620cb..000000000 --- a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/calculateStress.H +++ /dev/null @@ -1,53 +0,0 @@ - if (runTime.outputTime()) - { - // Displacement gradient - tetPointTensorField gradU = tetFec::grad(U); - - // Stress tensor - tetPointSymmTensorField sigma = - rho*(2.0*mu*symm(gradU) + lambda*I*tr(gradU)); - - // sigmaXX - tetPointScalarField sigmaXX - ( - IOobject - ( - "sigmaXX", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - sigma.component(symmTensor::XX) - ); - - // sigmaYY - tetPointScalarField sigmaYY - ( - IOobject - ( - "sigmaYY", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - sigma.component(symmTensor::YY) - ); - - // sigmaXY - tetPointScalarField sigmaXY - ( - IOobject - ( - "sigmaXY", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - sigma.component(symmTensor::XY) - ); - - runTime.write(); - } diff --git a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDEpsilonDSigma.H b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDEpsilonDSigma.H index dd03be9eb..50b2d2e95 100644 --- a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDEpsilonDSigma.H +++ b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDEpsilonDSigma.H @@ -1,5 +1,5 @@ -//- Increment of Green finite strain tensor +// Increment of Green finite strain tensor DEpsilon = symm(gradDU) + 0.5*symm((gradDU & gradU.T()) + (gradU & gradDU.T()) + (gradDU & gradDU.T())); -//- Increment of second Piola-Kirchhoff stress tensor +// Increment of second Piola-Kirchhoff stress tensor DSigma = 2*mu*(DEpsilon - DEpsilonP) + lambda*(I*tr(DEpsilon)); diff --git a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaExp.H b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaExp.H index c938e46ed..831730cf9 100644 --- a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaExp.H +++ b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaExp.H @@ -1,45 +1,46 @@ if(divDSigmaExpMethod == "standard") - { +{ divDSigmaExp = fvc::div - ( - (mu*gradDU.T()) - + (lambda*(I*tr(gradDU))) - - ((mu + lambda)*gradDU), - "div(sigma)" - ); - } - else if(divDSigmaExpMethod == "surface") - { - divDSigmaExp = - fvc::div( - mesh.magSf() - *( - muf*(n&fvc::interpolate(gradDU.T())) - + lambdaf*tr(fvc::interpolate(gradDU))*n - - (muf + lambdaf)*(n&fvc::interpolate(gradDU)) - ) - ); + ( + (mu*gradDU.T()) + + (lambda*(I*tr(gradDU))) + - ((mu + lambda)*gradDU), + "div(sigma)" + ); +} +else if(divDSigmaExpMethod == "surface") +{ + divDSigmaExp = fvc::div + ( + mesh.magSf()* + ( + muf*(n & fvc::interpolate(gradDU.T())) + + lambdaf*tr(fvc::interpolate(gradDU))*n + - (muf + lambdaf)*(n & fvc::interpolate(gradDU)) + ) + ); + // divDSigmaExp = fvc::div // ( // muf*(mesh.Sf() & fvc::interpolate(gradDU.T())) // + lambdaf*(mesh.Sf() & I*fvc::interpolate(tr(gradDU))) // - (muf + lambdaf)*(mesh.Sf() & fvc::interpolate(gradDU)) // ); - } - else if(divDSigmaExpMethod == "decompose") - { - surfaceTensorField shearGradDU = - ((I - n*n)&fvc::interpolate(gradDU)); +} +else if(divDSigmaExpMethod == "decompose") +{ + surfaceTensorField shearGradDU((I - n*n) & fvc::interpolate(gradDU)); - divDSigmaExp = - fvc::div( - mesh.magSf() - *( - - (muf + lambdaf)*(fvc::snGrad(DU)&(I - n*n)) - + lambdaf*tr(shearGradDU&(I - n*n))*n - + muf*(shearGradDU&n) - ) - ); + divDSigmaExp = fvc::div + ( + mesh.magSf()* + ( + - (muf + lambdaf)*(fvc::snGrad(DU) & (I - n*n)) + + lambdaf*tr(shearGradDU & (I - n*n))*n + + muf*(shearGradDU & n) + ) + ); + // divDSigmaExp = fvc::div // ( // mesh.magSf() @@ -49,19 +50,21 @@ if(divDSigmaExpMethod == "standard") // + muf*(shearGradDU&n) // ) // ); - } - else if(divDSigmaExpMethod == "laplacian") - { - divDSigmaExp = - - fvc::laplacian(mu + lambda, DU, "laplacian(DDU,DU)") - + fvc::div - ( - mu*gradDU.T() - + lambda*(I*tr(gradDU)), - "div(sigma)" - ); - } - else - { - FatalError << "divDSigmaExp method " << divDSigmaExpMethod << " not found!" << endl; - } +} +else if(divDSigmaExpMethod == "laplacian") +{ + divDSigmaExp = + - fvc::laplacian(mu + lambda, DU, "laplacian(DDU,DU)") + + fvc::div + ( + mu*gradDU.T() + + lambda*(I*tr(gradDU)), + "div(sigma)" + ); +} +else +{ + FatalErrorIn(args.executable()) + << "divDSigmaExp method " << divDSigmaExpMethod << " not found!" + << abort(FatalError); +} diff --git a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaNonLinExp.H b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaNonLinExp.H index 535509f4e..d95d0c730 100644 --- a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaNonLinExp.H +++ b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/calculateDivDSigmaNonLinExp.H @@ -1,32 +1,69 @@ if(divDSigmaNonLinExpMethod == "standard") - { +{ divDSigmaNonLinExp = fvc::div - ( - ( mu * ( - (gradDU & gradU.T()) - + (gradU & gradDU.T()) - + (gradDU & gradDU.T()) - ) ) - + ( lambda * 0.5 * tr( (gradDU & gradU.T()) - + (gradU & gradDU.T()) - + (gradDU & gradDU.T()) - ) * I ) - + ( DSigma & gradU ) - + ( (sigma + DSigma) & gradDU ), - "div(sigma)" - ); - } - else if(divDSigmaNonLinExpMethod == "surface") - { - divDSigmaNonLinExp = - fvc::div( - mesh.magSf() - *( - ( muf * (n & fvc::interpolate( (gradU & gradDU.T()) + (gradDU & gradU.T()) + (gradDU & gradDU.T()) )) ) - + ( 0.5*lambdaf * (n * tr(fvc::interpolate( (gradU & gradDU.T()) + (gradDU & gradU.T()) + (gradDU & gradDU.T()) ))) ) - + (n & fvc::interpolate( (DSigma & gradU) + ((sigma + DSigma) & gradDU) )) - ) - ); + ( + ( + mu* + ( + (gradDU & gradU.T()) + + (gradU & gradDU.T()) + + (gradDU & gradDU.T()) + ) + ) + + ( + 0.5*lambda* + tr + ( + (gradDU & gradU.T()) + + (gradU & gradDU.T()) + + (gradDU & gradDU.T()) + )*I + ) + + (DSigma & gradU) + + ((sigma + DSigma) & gradDU), + "div(sigma)" + ); +} +else if(divDSigmaNonLinExpMethod == "surface") +{ + divDSigmaNonLinExp = fvc::div + ( + mesh.magSf()* + ( + ( + muf* + ( + n & fvc::interpolate + ( + (gradU & gradDU.T()) + + (gradDU & gradU.T()) + + (gradDU & gradDU.T()) + ) + ) + ) + + ( + 0.5*lambdaf* + ( + n*tr + ( + fvc::interpolate + ( + (gradU & gradDU.T()) + + (gradDU & gradU.T()) + + (gradDU & gradDU.T()) + ) + ) + ) + ) + + ( + n & fvc::interpolate + ( + (DSigma & gradU) + + ((sigma + DSigma) & gradDU) + ) + ) + ) + ); // divDSigmaNonLinExp = // fvc::div @@ -48,8 +85,10 @@ if(divDSigmaNonLinExpMethod == "standard") // + ( (sigma + DSigma) & gradDU ) // ) ) // ); - } - else - { - FatalError << "divDSigmaExp method " << divDSigmaExpMethod << " not found!" << endl; - } +} +else +{ + FatalErrorIn(args.executable()) + << "divDSigmaExp method " << divDSigmaExpMethod << " not found!" + << abort(FatalError); +} diff --git a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/elasticPlasticNonLinTLSolidFoam.C b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/elasticPlasticNonLinTLSolidFoam.C index 1f2b203a9..5e24bf0be 100644 --- a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/elasticPlasticNonLinTLSolidFoam.C +++ b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/elasticPlasticNonLinTLSolidFoam.C @@ -66,114 +66,116 @@ int main(int argc, char *argv[]) scalar initialResidual = 0; lduMatrix::solverPerformance solverPerf; scalar relativeResidual = GREAT; - lduMatrix::debug=0; + lduMatrix::debug = 0; do - { - DU.storePrevIter(); + { + DU.storePrevIter(); # include "calculateDivDSigmaExp.H" # include "calculateDivDSigmaNonLinExp.H" - // Incremental form of the - // linear momentum conservation - // ensuring conservation of total momentum - fvVectorMatrix DUEqn - ( - fvm::d2dt2(rho, DU) - == - fvm::laplacian(2*muf + lambdaf, DU, "laplacian(DDU,DU)") - + divDSigmaExp - + divDSigmaNonLinExp - //- fvc::div(2*mu*DEpsilonP, "div(sigma)") - - fvc::div - ( - 2*muf*( mesh.Sf() & fvc::interpolate(DEpsilonP)) - ) - ); + // Incremental form of the + // linear momentum conservation + // ensuring conservation of total momentum + fvVectorMatrix DUEqn + ( + fvm::d2dt2(rho, DU) + == + fvm::laplacian(2*muf + lambdaf, DU, "laplacian(DDU,DU)") + + divDSigmaExp + + divDSigmaNonLinExp + //- fvc::div(2*mu*DEpsilonP, "div(sigma)") + - fvc::div + ( + 2*muf*(mesh.Sf() & fvc::interpolate(DEpsilonP)) + ) + ); - if (largeStrainOverRelax) - { - // the terms (gradDU & gradU.T()) and (gradU & gradDU.T()) - // are linearly dependent of DU and represent initial - // displacement effect - // which can cause convergence difficulties when treated - // explicitly - // so we implicitly over-relax with gradU & gradDU here - // which tends to help convergence - // this should improve convergence when gradU is large - // but maybe not execution time - DUEqn -= - fvm::laplacian - ( - (2*mu + lambda)*gradU, DU, "laplacian(DDU,DU)" - ) - - fvc::div( (2*mu + lambda)*(gradU&gradDU), "div(sigma)"); - } + if (largeStrainOverRelax) + { + // the terms (gradDU & gradU.T()) and (gradU & gradDU.T()) + // are linearly dependent of DU and represent initial + // displacement effect + // which can cause convergence difficulties when treated + // explicitly + // so we implicitly over-relax with gradU & gradDU here + // which tends to help convergence + // this should improve convergence when gradU is large + // but maybe not execution time + DUEqn -= + fvm::laplacian + ( + (2*mu + lambda)*gradU, DU, "laplacian(DDU,DU)" + ) + - fvc::div((2*mu + lambda)*(gradU & gradDU), "div(sigma)"); + } - if (nonLinearSemiImplicit) - { - // experimental - // we can treat the nonlinear term (gradDU & gradDU.T()) in a - // semi-implicit over-relaxed manner - // this should improve convergence when gradDU is large - // but maybe not execution time - DUEqn -= - fvm::laplacian - ( - (2*mu + lambda)*gradDU, DU, "laplacian(DDU,DU)" - ) - - fvc::div( (2*mu + lambda)*(gradDU&gradDU), "div(sigma)"); - } + if (nonLinearSemiImplicit) + { + // experimental + // we can treat the nonlinear term (gradDU & gradDU.T()) in a + // semi-implicit over-relaxed manner + // this should improve convergence when gradDU is large + // but maybe not execution time + DUEqn -= + fvm::laplacian + ( + (2*mu + lambda)*gradDU, DU, "laplacian(DDU,DU)" + ) + - fvc::div((2*mu + lambda)*(gradDU & gradDU), "div(sigma)"); + } - solverPerf = DUEqn.solve(); + solverPerf = DUEqn.solve(); - if (iCorr == 0) - { - initialResidual = solverPerf.initialResidual(); - } + if (iCorr == 0) + { + initialResidual = solverPerf.initialResidual(); + } - if (aitkenRelax) - { + if (aitkenRelax) + { # include "aitkenRelaxation.H" - } - else - { - DU.relax(); - } + } + else + { + DU.relax(); + } - gradDU = fvc::grad(DU); + gradDU = fvc::grad(DU); - // correct plasticty term - rheology.correct(); + // correct plasticty term + rheology.correct(); # include "calculateDEpsilonDSigma.H" # include "calculateRelativeResidual.H" - if (iCorr % infoFrequency == 0) - { - Info << "\tTime " << runTime.value() - << ", Corrector " << iCorr - << ", Solving for " << DU.name() - << " using " << solverPerf.solverName() - << ", res = " << solverPerf.initialResidual() - << ", rel res = " << relativeResidual; - if (aitkenRelax) - { - Info << ", aitken = " << aitkenTheta; - } - Info << ", iters = " << solverPerf.nIterations() << endl; - } - } + if (iCorr % infoFrequency == 0) + { + Info<< "\tTime " << runTime.value() + << ", Corrector " << iCorr + << ", Solving for " << DU.name() + << " using " << solverPerf.solverName() + << ", res = " << solverPerf.initialResidual() + << ", rel res = " << relativeResidual; + + if (aitkenRelax) + { + Info << ", aitken = " << aitkenTheta; + } + Info << ", iters = " << solverPerf.nIterations() << endl; + } + } while + ( + iCorr++ == 0 + || ( - iCorr++ == 0 - || - (//solverPerf.initialResidual() > convergenceTolerance - relativeResidual > convergenceTolerance - && - iCorr < nCorr) - ); + //solverPerf.initialResidual() > convergenceTolerance + relativeResidual > convergenceTolerance + && iCorr < nCorr + ) + ); Info<< nl << "Time " << runTime.value() << ", Solving for " << DU.name() << ", Initial residual = " << initialResidual diff --git a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/readDivDSigmaExpMethod.H b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/readDivDSigmaExpMethod.H index 1615c5211..421acb29d 100644 --- a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/readDivDSigmaExpMethod.H +++ b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/readDivDSigmaExpMethod.H @@ -1,9 +1,20 @@ //- how explicit component of sigma is to be calculated -word divDSigmaExpMethod(mesh.solutionDict().subDict("solidMechanics").lookup("divSigmaExp")); +word divDSigmaExpMethod +( + mesh.solutionDict().subDict("solidMechanics").lookup("divSigmaExp") +); + Info << "divSigmaExp method " << divDSigmaExpMethod << endl; -if(divDSigmaExpMethod != "standard" && divDSigmaExpMethod != "surface" && divDSigmaExpMethod != "decompose" && divDSigmaExpMethod != "laplacian") - { - FatalError << "divSigmaExp method " << divDSigmaExpMethod << " not found!" << nl - << "valid methods are:\nstandard\nsurface\ndecompose\nlaplacian" - << exit(FatalError); - } +if +( + divDSigmaExpMethod != "standard" + && divDSigmaExpMethod != "surface" + && divDSigmaExpMethod != "decompose" + && divDSigmaExpMethod != "laplacian" +) +{ + FatalErrorIn(args.executable()) + << "divSigmaExp method " << divDSigmaExpMethod << " not found!" << nl + << "valid methods are:\nstandard\nsurface\ndecompose\nlaplacian" + << exit(FatalError); +} diff --git a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/writeFields.H b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/writeFields.H index a3de4657b..c2db3454e 100644 --- a/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/writeFields.H +++ b/applications/solvers/solidMechanics/elasticPlasticNonLinTLSolidFoam/writeFields.H @@ -1,100 +1,95 @@ if (runTime.outputTime()) - { +{ volScalarField epsilonEq - ( - IOobject - ( - "epsilonEq", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - sqrt((2.0/3.0)*magSqr(dev(epsilon))) - ); - Info<< "Max epsilonEq = " << max(epsilonEq).value() - << endl; + ( + IOobject + ( + "epsilonEq", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + sqrt((2.0/3.0)*magSqr(dev(epsilon))) + ); + Info<< "Max epsilonEq = " << max(epsilonEq).value() << endl; volScalarField epsilonPEq - ( - IOobject - ( - "epsilonPEq", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE + ( + IOobject + ( + "epsilonPEq", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE ), - sqrt((2.0/3.0)*magSqr(dev(epsilonP))) - ); - Info<< "Max epsilonPEq = " << max(epsilonPEq).value() - << endl; + sqrt((2.0/3.0)*magSqr(dev(epsilonP))) + ); + Info<< "Max epsilonPEq = " << max(epsilonPEq).value()<< endl; volScalarField sigmaEq - ( - IOobject - ( - "sigmaEq", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - sqrt((3.0/2.0)*magSqr(dev(sigma))) - ); - - Info<< "Max sigmaEq = " << max(sigmaEq).value() - << endl; + ( + IOobject + ( + "sigmaEq", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + sqrt((3.0/2.0)*magSqr(dev(sigma))) + ); - // deformation gradient + Info<< "Max sigmaEq = " << max(sigmaEq).value() << endl; + + // Deformation gradient volTensorField F = I + gradU; volScalarField J = det(F); - //- Calculate Cauchy stress + // Calculate Cauchy stress volSymmTensorField sigmaCauchy - ( - IOobject - ( - "sigmaCauchy", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - (1/J) * symm(F.T() & sigma & F) - ); + ( + IOobject + ( + "sigmaCauchy", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + (1/J)*symm(F.T() & sigma & F) + ); - //- Cauchy von Mises stress + // Cauchy von Mises stress volScalarField sigmaCauchyEq - ( - IOobject - ( - "sigmaCauchyEq", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::AUTO_WRITE - ), - sqrt((3.0/2.0)*magSqr(dev(sigmaCauchy))) - ); - - Info<< "Max sigmaCauchyEq = " << max(sigmaCauchyEq).value() - << endl; + ( + IOobject + ( + "sigmaCauchyEq", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + sqrt((3.0/2.0)*magSqr(dev(sigmaCauchy))) + ); - //volTensorField Finv = inv(F); - // volSymmTensorField epsilonAlmansi -// ( -// IOobject -// ( -// "epsilonAlmansi", -// runTime.timeName(), -// mesh, -// IOobject::NO_READ, -// IOobject::AUTO_WRITE -// ), -// symm(Finv & epsilon & Finv.T()) -// ); + Info<< "Max sigmaCauchyEq = " << max(sigmaCauchyEq).value() << endl; +// volTensorField Finv = inv(F); +// volSymmTensorField epsilonAlmansi +// ( +// IOobject +// ( +// "epsilonAlmansi", +// runTime.timeName(), +// mesh, +// IOobject::NO_READ, +// IOobject::AUTO_WRITE +// ), +// symm(Finv & epsilon & Finv.T()) +// ); runTime.write(); - } +} diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/Make/files b/applications/solvers/solidMechanics/stressFemFoam/Make/files similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/Make/files rename to applications/solvers/solidMechanics/stressFemFoam/Make/files diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/Make/options b/applications/solvers/solidMechanics/stressFemFoam/Make/options similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/Make/options rename to applications/solvers/solidMechanics/stressFemFoam/Make/options diff --git a/applications/solvers/solidMechanics/stressFemFoam/calculateStress.H b/applications/solvers/solidMechanics/stressFemFoam/calculateStress.H new file mode 100644 index 000000000..3e4619cef --- /dev/null +++ b/applications/solvers/solidMechanics/stressFemFoam/calculateStress.H @@ -0,0 +1,120 @@ + if (runTime.outputTime()) + { + // Displacement gradient + tetPointTensorField gradU = tetFec::grad(U); + + // Stress tensor + tetPointSymmTensorField sigma = + rho*(2.0*mu*symm(gradU) + lambda*I*tr(gradU)); + + + // Create pointMesh for field post-processing + const pointMesh& pMesh = pointMesh::New(mesh); + + // U + pointVectorField Up + ( + IOobject + ( + "Up", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + pMesh, + U.dimensions() + ); + + Up.internalField() = vectorField::subField + ( + U.internalField(), + pMesh.size() + ); + + // sigmaEq + pointScalarField sigmaEq + ( + IOobject + ( + "sigmaEq", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + pMesh, + sigma.dimensions() + ); + + sigmaEq.internalField() = scalarField::subField + ( + sqrt((3.0/2.0)*magSqr(dev(sigma.internalField())))(), + pMesh.size() + ); + + // sigmaXX + pointScalarField sigmaXX + ( + IOobject + ( + "sigmaXX", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + pMesh, + sigma.dimensions() + ); + + sigmaXX.internalField() = scalarField::subField + ( + sigma.component(symmTensor::XX)().internalField(), + pMesh.size() + ); + + // sigmaYY + pointScalarField sigmaYY + ( + IOobject + ( + "sigmaYY", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + pMesh, + sigma.dimensions() + ); + + sigmaYY.internalField() = scalarField::subField + ( + sigma.component(symmTensor::YY)().internalField(), + pMesh.size() + ); + + // sigmaXY + pointScalarField sigmaXY + ( + IOobject + ( + "sigmaXY", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::AUTO_WRITE + ), + pMesh, + sigma.dimensions() + ); + + sigmaXY.internalField() = scalarField::subField + ( + sigma.component(symmTensor::XY)().internalField(), + pMesh.size() + ); + + runTime.write(); + } diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/createFemFields.H b/applications/solvers/solidMechanics/stressFemFoam/createFemFields.H similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/createFemFields.H rename to applications/solvers/solidMechanics/stressFemFoam/createFemFields.H diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/createFemMesh.H b/applications/solvers/solidMechanics/stressFemFoam/createFemMesh.H similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/createFemMesh.H rename to applications/solvers/solidMechanics/stressFemFoam/createFemMesh.H diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/createFields.H b/applications/solvers/solidMechanics/stressFemFoam/createFields.H similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/createFields.H rename to applications/solvers/solidMechanics/stressFemFoam/createFields.H diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.C b/applications/solvers/solidMechanics/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.C similarity index 96% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.C rename to applications/solvers/solidMechanics/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.C index b5cf8cb66..a10dadc9d 100644 --- a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.C +++ b/applications/solvers/solidMechanics/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.C @@ -46,8 +46,9 @@ template class PointPatch, template class MatrixType > -void TractionPointPatchVectorField -::checkFieldSize() const +void +TractionPointPatchVectorField:: +checkFieldSize() const { if ( @@ -229,7 +230,7 @@ addBoundarySourceDiag const vectorField& pointNormals = this->patch().pointNormals(); const label nFaces = this->patch().nFaces(); - vectorField& source = matrix.source(); + vectorField& b = matrix.b(); for (label faceI = 0; faceI < nFaces; faceI++) { @@ -240,7 +241,7 @@ addBoundarySourceDiag { const triFace& tri = faceTriangles[triI]; - source[meshPoints[tri[0]]] += + b[meshPoints[tri[0]]] += tri.mag(points)* ( traction_[faceI]/6.0 @@ -251,7 +252,7 @@ addBoundarySourceDiag - pressure_[faceI]*pointNormals[tri[2]]/12.0 ); - source[meshPoints[tri[1]]] += + b[meshPoints[tri[1]]] += tri.mag(points)* ( traction_[faceI]/6.0 @@ -262,7 +263,7 @@ addBoundarySourceDiag - pressure_[faceI]*pointNormals[tri[0]]/12.0 ); - source[meshPoints[tri[2]]] += + b[meshPoints[tri[2]]] += tri.mag(points)* ( traction_[faceI]/6.0 diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.H b/applications/solvers/solidMechanics/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.H similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.H rename to applications/solvers/solidMechanics/stressFemFoam/femStress/Traction/TractionPointPatchVectorField.H diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/tractionTetPolyPatchVectorField.C b/applications/solvers/solidMechanics/stressFemFoam/femStress/tractionTetPolyPatchVectorField.C similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/tractionTetPolyPatchVectorField.C rename to applications/solvers/solidMechanics/stressFemFoam/femStress/tractionTetPolyPatchVectorField.C diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/tractionTetPolyPatchVectorField.H b/applications/solvers/solidMechanics/stressFemFoam/femStress/tractionTetPolyPatchVectorField.H similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/femStress/tractionTetPolyPatchVectorField.H rename to applications/solvers/solidMechanics/stressFemFoam/femStress/tractionTetPolyPatchVectorField.H diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/readControls.H b/applications/solvers/solidMechanics/stressFemFoam/readControls.H similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/readControls.H rename to applications/solvers/solidMechanics/stressFemFoam/readControls.H diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/readMechanicalProperties.H b/applications/solvers/solidMechanics/stressFemFoam/readMechanicalProperties.H similarity index 100% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/readMechanicalProperties.H rename to applications/solvers/solidMechanics/stressFemFoam/readMechanicalProperties.H diff --git a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/stressFemFoam.C b/applications/solvers/solidMechanics/stressFemFoam/stressFemFoam.C similarity index 97% rename from applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/stressFemFoam.C rename to applications/solvers/solidMechanics/stressFemFoam/stressFemFoam.C index 65b8fae2a..50daab29a 100644 --- a/applications/solvers/solidMechanics/deprecatedSolvers/stressFemFoam/stressFemFoam.C +++ b/applications/solvers/solidMechanics/stressFemFoam/stressFemFoam.C @@ -36,6 +36,8 @@ Description #include "fvCFD.H" +#include "pointMesh.H" +#include "pointFields.H" #include "tetPolyMesh.H" #include "tetPointFields.H" #include "tetFem.H" @@ -69,7 +71,7 @@ int main(int argc, char *argv[]) + tetFem::laplacianTrace(lambda, U) ); - solve(UEqn); + UEqn.solve(); # include "calculateStress.H" diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/Allclean b/tutorials/solidMechanics/stressFemFoam/Allclean similarity index 100% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/Allclean rename to tutorials/solidMechanics/stressFemFoam/Allclean diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/Allrun b/tutorials/solidMechanics/stressFemFoam/Allrun similarity index 100% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/Allrun rename to tutorials/solidMechanics/stressFemFoam/Allrun diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/0/U b/tutorials/solidMechanics/stressFemFoam/plateHole/0/U similarity index 100% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/0/U rename to tutorials/solidMechanics/stressFemFoam/plateHole/0/U diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/mechanicalProperties b/tutorials/solidMechanics/stressFemFoam/plateHole/constant/mechanicalProperties similarity index 100% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/mechanicalProperties rename to tutorials/solidMechanics/stressFemFoam/plateHole/constant/mechanicalProperties diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/polyMesh/blockMeshDict b/tutorials/solidMechanics/stressFemFoam/plateHole/constant/polyMesh/blockMeshDict similarity index 72% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/polyMesh/blockMeshDict rename to tutorials/solidMechanics/stressFemFoam/plateHole/constant/polyMesh/blockMeshDict index 8b27dc96d..227250905 100644 --- a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/polyMesh/blockMeshDict +++ b/tutorials/solidMechanics/stressFemFoam/plateHole/constant/polyMesh/blockMeshDict @@ -19,28 +19,28 @@ convertToMeters 1; vertices ( - (0.5 0 0) - (1 0 0) - (2 0 0) - (2 0.707107 0) - (0.707107 0.707107 0) - (0.353553 0.353553 0) - (2 2 0) - (0.707107 2 0) - (0 2 0) - (0 1 0) - (0 0.5 0) - (0.5 0 0.5) - (1 0 0.5) - (2 0 0.5) - (2 0.707107 0.5) - (0.707107 0.707107 0.5) - (0.353553 0.353553 0.5) - (2 2 0.5) - (0.707107 2 0.5) - (0 2 0.5) - (0 1 0.5) - (0 0.5 0.5) + (0.5 0 -0.1) + (1 0 -0.1) + (2 0 -0.1) + (2 0.707107 -0.1) + (0.707107 0.707107 -0.1) + (0.353553 0.353553 -0.1) + (2 2 -0.1) + (0.707107 2 -0.1) + (0 2 -0.1) + (0 1 -0.1) + (0 0.5 -0.1) + (0.5 0 0.1) + (1 0 0.1) + (2 0 0.1) + (2 0.707107 0.1) + (0.707107 0.707107 0.1) + (0.353553 0.353553 0.1) + (2 2 0.1) + (0.707107 2 0.1) + (0 2 0.1) + (0 1 0.1) + (0 0.5 0.1) ); blocks @@ -54,14 +54,14 @@ blocks edges ( - arc 0 5 (0.469846 0.17101 0) - arc 5 10 (0.17101 0.469846 0) - arc 1 4 (0.939693 0.34202 0) - arc 4 9 (0.34202 0.939693 0) - arc 11 16 (0.469846 0.17101 0.5) - arc 16 21 (0.17101 0.469846 0.5) - arc 12 15 (0.939693 0.34202 0.5) - arc 15 20 (0.34202 0.939693 0.5) + arc 0 5 (0.469846 0.17101 -0.1) + arc 5 10 (0.17101 0.469846 -0.1) + arc 1 4 (0.939693 0.34202 -0.1) + arc 4 9 (0.34202 0.939693 -0.1) + arc 11 16 (0.469846 0.17101 0.1) + arc 16 21 (0.17101 0.469846 0.1) + arc 12 15 (0.939693 0.34202 0.1) + arc 15 20 (0.34202 0.939693 0.1) ); patches diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/polyMesh/boundary b/tutorials/solidMechanics/stressFemFoam/plateHole/constant/polyMesh/boundary similarity index 96% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/polyMesh/boundary rename to tutorials/solidMechanics/stressFemFoam/plateHole/constant/polyMesh/boundary index f7b45bdb7..1401a1b37 100644 --- a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/constant/polyMesh/boundary +++ b/tutorials/solidMechanics/stressFemFoam/plateHole/constant/polyMesh/boundary @@ -1,7 +1,7 @@ /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | foam-extend: Open Source CFD | -| \\ / O peration | Version: 3.0 | +| \\ / O peration | Version: 3.0 | | \\ / A nd | Web: http://www.extend-project.de | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/controlDict b/tutorials/solidMechanics/stressFemFoam/plateHole/system/controlDict similarity index 94% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/controlDict rename to tutorials/solidMechanics/stressFemFoam/plateHole/system/controlDict index 8a9087959..80f2e4a5f 100644 --- a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/controlDict +++ b/tutorials/solidMechanics/stressFemFoam/plateHole/system/controlDict @@ -15,7 +15,7 @@ FoamFile } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -application stressFemFoam; +application stressFemFoam; startFrom startTime; @@ -23,13 +23,13 @@ startTime 0; stopAt endTime; -endTime 100; +endTime 1; deltaT 1; writeControl timeStep; -writeInterval 20; +writeInterval 1; purgeWrite 0; diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/fvSchemes b/tutorials/solidMechanics/stressFemFoam/plateHole/system/fvSchemes similarity index 100% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/fvSchemes rename to tutorials/solidMechanics/stressFemFoam/plateHole/system/fvSchemes diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/fvSolution b/tutorials/solidMechanics/stressFemFoam/plateHole/system/fvSolution similarity index 100% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/fvSolution rename to tutorials/solidMechanics/stressFemFoam/plateHole/system/fvSolution diff --git a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/tetFemSolution b/tutorials/solidMechanics/stressFemFoam/plateHole/system/tetFemSolution similarity index 86% rename from tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/tetFemSolution rename to tutorials/solidMechanics/stressFemFoam/plateHole/system/tetFemSolution index 91a9a7218..85927aea3 100644 --- a/tutorials/solidMechanics/deprecatedTutorials/stressFemFoam/plateHole/system/tetFemSolution +++ b/tutorials/solidMechanics/stressFemFoam/plateHole/system/tetFemSolution @@ -19,11 +19,15 @@ solvers { U { - solver PCG; - preconditioner DIC; + solver CG; + preconditioner Cholesky; + minIter 1; + maxIter 2000; + tolerance 1e-06; - relTol 0.1; + relTol 0; } } + // ************************************************************************* //