diff --git a/applications/solvers/solidMechanics/elasticThermalSolidFoam/createFields.H b/applications/solvers/solidMechanics/elasticThermalSolidFoam/createFields.H index e789e67f2..0ab66ab7e 100644 --- a/applications/solvers/solidMechanics/elasticThermalSolidFoam/createFields.H +++ b/applications/solvers/solidMechanics/elasticThermalSolidFoam/createFields.H @@ -12,21 +12,20 @@ mesh ); - // volTensorField gradU = fvc::grad(U); volTensorField gradU ( IOobject ( - "grad(U)", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::NO_WRITE - ), + "grad(U)", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), mesh, dimensionedTensor("zero", dimless, tensor::zero) - ); - //surfaceVectorField snGradU = fvc::snGrad(U); + ); + surfaceVectorField snGradU ( IOobject @@ -38,7 +37,7 @@ IOobject::NO_WRITE ), mesh, - dimensionedVector("zero", dimless, vector::zero) + dimensionedVector("zero", dimless, vector::zero) ); volSymmTensorField epsilon @@ -115,31 +114,33 @@ thermalModel thermal(T); volScalarField C = thermal.C(); - volScalarField k = thermal.k(); + volScalarField k + ( + "DT", + thermal.k() + ); volScalarField threeKalpha = rheology.threeK()*rho*thermal.alpha(); - surfaceScalarField threeKalphaf = fvc::interpolate(threeKalpha, "threeKalpha"); + surfaceScalarField threeKalphaf = + fvc::interpolate(threeKalpha, "threeKalpha"); volScalarField T0 = thermal.T0(); volScalarField rhoC = rho*C; -// for aitken relaxation -volVectorField aitkenDelta - ( + // for aitken relaxation + volVectorField aitkenDelta + ( IOobject ( - "aitkenDelta", - runTime.timeName(), - mesh, - IOobject::NO_READ, - IOobject::NO_WRITE - ), + "aitkenDelta", + runTime.timeName(), + mesh, + IOobject::NO_READ, + IOobject::NO_WRITE + ), mesh, dimensionedVector("zero", dimLength, vector::zero) ); -// aitken relaxation factor -scalar aitkenInitialRes = 1.0; -scalar aitkenTheta = 0.01; -// if(mesh.relax(U.name())) -// { -// aitkenTheta = mesh.relaxationFactor(U.name()); -// } + + // aitken relaxation factor + scalar aitkenInitialRes = 1.0; + scalar aitkenTheta = 0.01; diff --git a/applications/solvers/solidMechanics/elasticThermalSolidFoam/elasticThermalSolidFoam.C b/applications/solvers/solidMechanics/elasticThermalSolidFoam/elasticThermalSolidFoam.C index 3fe868901..c90b1096e 100644 --- a/applications/solvers/solidMechanics/elasticThermalSolidFoam/elasticThermalSolidFoam.C +++ b/applications/solvers/solidMechanics/elasticThermalSolidFoam/elasticThermalSolidFoam.C @@ -46,100 +46,99 @@ Author int main(int argc, char *argv[]) { -# include "setRootCase.H" -# include "createTime.H" -# include "createMesh.H" -# include "createFields.H" -# include "readDivSigmaExpMethod.H" +# include "setRootCase.H" +# include "createTime.H" +# include "createMesh.H" +# include "createFields.H" +# include "readDivSigmaExpMethod.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // - Info<< "\nStarting time loop\n" << endl; + Info<< "\nStarting time loop\n" << endl; - while(runTime.loop()) + while(runTime.loop()) { - Info<< "Time: " << runTime.timeName() << nl << endl; + Info<< "Time: " << runTime.timeName() << nl << endl; -# include "readSolidMechanicsControls.H" +# include "readSolidMechanicsControls.H" - int iCorr = 0; - scalar initialResidual = 1.0; - scalar relResT = 1.0; - scalar relResU = 1.0; - lduMatrix::solverPerformance solverPerfU; - lduMatrix::solverPerformance solverPerfT; - lduMatrix::debug = 0; + int iCorr = 0; + scalar initialResidual = 1.0; + scalar relResT = 1.0; + scalar relResU = 1.0; + lduMatrix::solverPerformance solverPerfU; + lduMatrix::solverPerformance solverPerfT; + lduMatrix::debug = 0; - // solve energy equation for temperature - // the loop is for non-orthogonal corrections - Info<< "Solving for " << T.name() << nl; - do - { - T.storePrevIter(); + // solve energy equation for temperature + // the loop is for non-orthogonal corrections + Info<< "Solving for " << T.name() << nl; + do + { + T.storePrevIter(); - fvScalarMatrix TEqn - ( - rhoC*fvm::ddt(T) == fvm::laplacian(k, T, "laplacian(k,T)") - ); + fvScalarMatrix TEqn + ( + rhoC*fvm::ddt(T) == fvm::laplacian(k, T, "laplacian(k,T)") + ); - solverPerfT = TEqn.solve(); + solverPerfT = TEqn.solve(); - T.relax(); + T.relax(); -# include "calculateRelResT.H" +# include "calculateRelResT.H" - if (iCorr % infoFrequency == 0) - { - Info<< "\tCorrector " << iCorr - << ", residual = " << solverPerfT.initialResidual() - << ", relative res = " << relResT - << ", inner iters = " << solverPerfT.nIterations() << endl; - } - } - while - ( - relResT > convergenceToleranceT - && - ++iCorr < nCorr - ); + if (iCorr % infoFrequency == 0) + { + Info<< "\tCorrector " << iCorr + << ", residual = " << solverPerfT.initialResidual() + << ", relative res = " << relResT + << ", inner iters = " << solverPerfT.nIterations() << endl; + } + } + while + ( + relResT > convergenceToleranceT + && ++iCorr < nCorr + ); - Info<< "Solved for " << T.name() - << " using " << solverPerfT.solverName() - << " in " << iCorr << " iterations" - << ", residual = " << solverPerfT.initialResidual() - << ", relative res = " << relResT << nl - << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" - << ", ClockTime = " << runTime.elapsedClockTime() << " s" - << endl; + Info<< "Solved for " << T.name() + << " using " << solverPerfT.solverName() + << " in " << iCorr << " iterations" + << ", residual = " << solverPerfT.initialResidual() + << ", relative res = " << relResT << nl + << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << ", ClockTime = " << runTime.elapsedClockTime() << " s" + << endl; - // solve momentum equation for displacement - iCorr = 0; - volVectorField gradThreeKalphaDeltaT = - fvc::grad(threeKalpha*(T-T0), "grad(threeKalphaDeltaT)"); - surfaceVectorField threeKalphaDeltaTf = - mesh.Sf()*threeKalphaf*fvc::interpolate(T-T0, "deltaT"); + // Solve momentum equation for displacement + iCorr = 0; + volVectorField gradThreeKalphaDeltaT = + fvc::grad(threeKalpha*(T-T0), "grad(threeKalphaDeltaT)"); + surfaceVectorField threeKalphaDeltaTf = + mesh.Sf()*threeKalphaf*fvc::interpolate(T-T0, "deltaT"); - Info<< "Solving for " << U.name() << nl; - do + Info<< "Solving for " << U.name() << nl; + do { U.storePrevIter(); -# include "calculateDivSigmaExp.H" +# include "calculateDivSigmaExp.H" // Linear momentum equaiton fvVectorMatrix UEqn - ( - rho*fvm::d2dt2(U) - == - fvm::laplacian(2*muf + lambdaf, U, "laplacian(DU,U)") - + divSigmaExp - ); + ( + rho*fvm::d2dt2(U) + == + fvm::laplacian(2*muf + lambdaf, U, "laplacian(DU,U)") + + divSigmaExp + ); solverPerfU = UEqn.solve(); if (aitkenRelax) { -# include "aitkenRelaxation.H" +# include "aitkenRelaxation.H" } else { @@ -148,7 +147,7 @@ int main(int argc, char *argv[]) gradU = fvc::grad(U); -# include "calculateRelResU.H" +# include "calculateRelResU.H" if (iCorr == 0) { @@ -160,6 +159,7 @@ int main(int argc, char *argv[]) Info<< "\tCorrector " << iCorr << ", residual = " << solverPerfU.initialResidual() << ", relative res = " << relResU; + if (aitkenRelax) { Info << ", aitken = " << aitkenTheta; @@ -167,25 +167,24 @@ int main(int argc, char *argv[]) Info<< ", inner iters = " << solverPerfU.nIterations() << endl; } } - while - ( - iCorr++ == 0 - || - (//solverPerfU.initialResidual() > convergenceTolerance - relResU > convergenceToleranceU - && - iCorr < nCorr) - ); + while + ( + iCorr++ == 0 + || ( + relResU > convergenceToleranceU + && iCorr < nCorr + ) + ); - Info<< "Solved for " << U.name() - << " using " << solverPerfU.solverName() - << " in " << iCorr << " iterations" - << ", initial res = " << initialResidual - << ", final res = " << solverPerfU.initialResidual() - << ", final rel res = " << relResU << nl - << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" - << ", ClockTime = " << runTime.elapsedClockTime() << " s" - << endl; + Info<< "Solved for " << U.name() + << " using " << solverPerfU.solverName() + << " in " << iCorr << " iterations" + << ", initial res = " << initialResidual + << ", final res = " << solverPerfU.initialResidual() + << ", final rel res = " << relResU << nl + << "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << ", ClockTime = " << runTime.elapsedClockTime() << " s" + << endl; # include "calculateEpsilonSigma.H" # include "writeFields.H" diff --git a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/T b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/T index a693bb903..4824e3cc0 100644 --- a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/T +++ b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/T @@ -32,27 +32,23 @@ boundaryField inside { - type fixedValue; - value uniform 100; + type fixedValue; + value uniform 100; } - outside { - type fixedValue; - value uniform 100; + type fixedValue; + value uniform 100; } - front { - type empty; - //type symmetryPlane; + type empty; } - back { - type empty; - //type symmetryPlane; + type empty; } } + // ************************************************************************* // diff --git a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/U b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/U index cbe770106..dd1957b17 100644 --- a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/U +++ b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/0/U @@ -29,32 +29,31 @@ boundaryField { type symmetryPlane; } - inside { - type solidTraction; - traction uniform ( 0 0 0 ); - pressure uniform 50e6; - value uniform (0 0 0); + type solidTraction; + traction uniform ( 0 0 0 ); + pressure uniform 50e6; + DT k; + value uniform (0 0 0); } - outside { type solidTraction; traction uniform ( 0 0 0 ); pressure uniform 0.1e6; + DT k; value uniform (0 0 0); } - front { - type empty; + type empty; } - back { - type empty; + type empty; } } + // ************************************************************************* // diff --git a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/Allrun b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/Allrun index 1eb6aec58..394b0dcf8 100755 --- a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/Allrun +++ b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/Allrun @@ -9,4 +9,4 @@ runApplication blockMesh runApplication $application (cd analyticalHotCylinder && runApplication wmake) -runApplication analyticalHotCylinder \ No newline at end of file +runApplication analyticalHotCylinder diff --git a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/polyMesh/blockMeshDict b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/polyMesh/blockMeshDict index 5cd413a66..e0274395f 100644 --- a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/polyMesh/blockMeshDict +++ b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/polyMesh/blockMeshDict @@ -19,15 +19,14 @@ convertToMeters 1; vertices ( - (0 0.5 0) - (0.5 0 0) - (0.7 0 0) - (0 0.7 0) - (0 0.5 0.1) - (0.5 0 0.1) - (0.7 0 0.1) - (0 0.7 0.1) - + (0 0.5 0) + (0.5 0 0) + (0.7 0 0) + (0 0.7 0) + (0 0.5 0.1) + (0.5 0 0.1) + (0.7 0 0.1) + (0 0.7 0.1) ); blocks @@ -45,34 +44,30 @@ edges patches ( - + symmetryPlane left + ( + (4 7 3 0) + ) + symmetryPlane bottom + ( + (1 2 6 5) + ) + patch inside + ( + (0 1 5 4) + ) + patch outside + ( + (7 6 2 3) + ) empty back ( - (3 2 1 0) + (3 2 1 0) ) empty front ( - (4 5 6 7) - ) - - symmetryPlane left - ( - (4 7 3 0) - ) - - symmetryPlane bottom - ( - (1 2 6 5) - ) - - patch inside - ( - (0 1 5 4) - ) - patch outside - ( - (7 6 2 3) - ) + (4 5 6 7) + ) ); mergePatchPairs diff --git a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/polyMesh/boundary b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/polyMesh/boundary index 498479cec..44af09850 100644 --- a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/polyMesh/boundary +++ b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/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 | | \*---------------------------------------------------------------------------*/ @@ -17,41 +17,41 @@ FoamFile 6 ( - back - { - type empty; - nFaces 600; - startFace 1130; - } - front - { - type empty; - nFaces 600; - startFace 1730; - } left { type symmetryPlane; nFaces 10; - startFace 2330; + startFace 1130; } bottom { type symmetryPlane; nFaces 10; - startFace 2340; + startFace 1140; } inside { type patch; nFaces 60; - startFace 2350; + startFace 1150; } outside { type patch; nFaces 60; - startFace 2410; + startFace 1210; + } + back + { + type empty; + nFaces 600; + startFace 1270; + } + front + { + type empty; + nFaces 600; + startFace 1870; } ) diff --git a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/thermalProperties b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/thermalProperties index 41a212479..a9d516ba9 100644 --- a/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/thermalProperties +++ b/tutorials/solidMechanics/elasticThermalSolidFoam/hotCylinder/constant/thermalProperties @@ -16,11 +16,12 @@ FoamFile thermal { - type constant; - C C [0 2 -2 -1 0 0 0] 434; - k k [1 1 -3 -1 0 0 0] 250; - alpha alpha [0 0 0 -1 0 0 0] 2.3e-05; - T0 T0 [0 0 0 1 0 0 0] 0; + type constant; + C C [0 2 -2 -1 0 0 0] 434; + k k [1 1 -3 -1 0 0 0] 250; + alpha alpha [0 0 0 -1 0 0 0] 2.3e-05; + T0 T0 [0 0 0 1 0 0 0] 0; } + // ************************************************************************* //