diff --git a/.gitignore b/.gitignore index f567511db3..3cc18eb468 100644 --- a/.gitignore +++ b/.gitignore @@ -55,6 +55,10 @@ doc/[Dd]oxygen/man /*.html /doc/*.html +# untracked configuration files +/etc/prefs.csh +/etc/prefs.sh + # source packages - anywhere *.tar.bz2 *.tar.gz diff --git a/README b/README index df56c9c9a4..dbb820aad3 100644 --- a/README +++ b/README @@ -2,7 +2,7 @@ # #+TITLE: OpenFOAM README for version 1.6 #+AUTHOR: OpenCFD Ltd. -#+DATE: November 2009 +#+DATE: March 2010 #+LINK: http://www.opencfd.co.uk #+OPTIONS: author:nil ^:{} @@ -23,10 +23,10 @@ section "Running OpenFOAM in 32-bit mode". *** Qt (from http://trolltech.com/products/qt) - The ParaView 3.6.1 visualisation package requires Qt to be installed on the + The ParaView 3.7.0 visualisation package requires Qt to be installed on the system. ParaView's producers state that ParaView is only officially - supported on Qt version 4.3.x. However, we have found in limited tests that - ParaView works satisfactorily with newer versions of Qt than 4.3.x. To + supported on Qt version 4.6.x. However, we have found in limited tests that + ParaView works satisfactorily with Qt than 4.5.x. To check whether Qt4 is installed, and the version, type: + qmake --version @@ -44,13 +44,14 @@ + openSUSE-10.3: Version 4.3.1 + openSUSE-11.0: Version 4.4.0 + openSUSE-11.1: Version 4.4.3 + + openSUSE-11.2: Version 4.5.3 Compilation and running of ParaView has been successful using the libraries downloaded in the "libqt4-dev" package on ubuntu. If you don't have an appropriate version of Qt installed you can download - the sources from TrollTech e.g.: - ftp://ftp.trolltech.com/qt/source/qt-x11-opensource-src-4.3.5.tar.bz2 + the sources e.g.: + http://get.qt.nokia.com/qt/source/qt-everywhere-opensource-src-4.6.2.tar.gz and compile and install in /usr/local or some other location that does not conflict with the pre-installed version. @@ -108,11 +109,11 @@ * Building from Sources (Optional) If you cannot find an appropriate binary pack for your platform, you can build the complete OpenFOAM from the source-pack. You will first need to compile or - obtain a recent version of gcc (we recomend gcc-4.3.?) for your platform, + obtain a recent version of gcc (we recommend gcc-4.4.?) for your platform, which may be obtained from http://gcc.gnu.org/. Install the compiler in - $WM_THIRD_PARTY_DIR/gcc-/platforms/$WM_ARCH$WM_COMPILER_ARCH/ + $WM_THIRD_PARTY_DIR/platforms/$WM_ARCH$WM_COMPILER_ARCH/gcc- and change the gcc version number in $WM_PROJECT_DIR/etc/settings.sh and $WM_PROJECT_DIR/etc/settings.csh appropriately and finally update the environment variables as in section 3. @@ -157,15 +158,16 @@ Refer to the OpenFOAM User Guide at http://www.OpenFOAM.org/doc/user.html for more information. -* Compiling Paraview 3.6.1 and the PV3FoamReader module +* Compiling Paraview 3.7.0 and the PV3FoamReader module If there are problems encountered with ParaView, then it may be necessary to compile ParaView from sources. The compilation is a fairly simple process using the makeParaView script (found in ThirdParty directory), which has worked in our tests with other - packages supplied in the ThirdParty directory, namely cmake-2.6.4 and - gcc-4.3.3. Execute the following: + packages supplied in the ThirdParty directory, namely cmake-2.8.0 and + gcc-4.4.3. Execute the following: + cd $WM_THIRD_PARTY_DIR - + rm -rf paraview-3.6.1/platforms + + rm -rf paraview-3.7.0/platforms + + rm -rf platforms/*/paraview-3.7.0 + ./makeParaView The PV3blockMeshReader and the PV3FoamReader ParaView plugins are compiled @@ -177,7 +179,7 @@ *** Compiling Paraview with a local version of Qt If the user still encounters problems with ParaView, it may relate to the version of Qt, in which case, it is recommended that the user first - downloads a supported version of Qt /e.g./ 4.3.5 as described in the section + downloads a supported version of Qt /e.g./ 4.5.3 as described in the section on "Qt". The user should unpack the source pack in the $WM_THIRD_PARTY_DIR. Then the user can build Qt by executing from within $WM_THIRD_PARTY_DIR: + ./makeQt diff --git a/TODO b/TODO index c2fd049c26..5c801b5377 100644 --- a/TODO +++ b/TODO @@ -24,6 +24,8 @@ OK - amg. Tested on unitTestCases/singleCyclic/ - initTransfer in GAMGprocessorInterfaces using nonblocking+tags +untested. + - test createPatch pointSync - pointFields on cyclics. volPointInterpolation. - jumpCyclics diff --git a/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C b/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C index aa9fd58b70..fdddddfbfe 100644 --- a/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C +++ b/applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.C @@ -270,7 +270,7 @@ void PDRkEpsilon::correct() } tmp tgradU = fvc::grad(U_); - volScalarField G = 2*mut_*(tgradU() && dev(symm(tgradU()))); + volScalarField G("RASModel::G", mut_*(tgradU() && dev(twoSymm(tgradU())))); tgradU.clear(); // Update espsilon and G at the wall diff --git a/applications/solvers/combustion/fireFoam/combustionModels/Make/options b/applications/solvers/combustion/fireFoam/combustionModels/Make/options index 42f59c8e17..3f557113c4 100644 --- a/applications/solvers/combustion/fireFoam/combustionModels/Make/options +++ b/applications/solvers/combustion/fireFoam/combustionModels/Make/options @@ -1,6 +1,4 @@ EXE_INC = \ - -I../sensibleEnthalpyCombustionThermophysicalModels/basic/lnInclude \ - -I../sensibleEnthalpyCombustionThermophysicalModels/reactionThermo/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \ -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \ diff --git a/applications/solvers/compressible/rhoPorousMRFPimpleFoam/Make/files b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/Make/files new file mode 100644 index 0000000000..a798455eb7 --- /dev/null +++ b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/Make/files @@ -0,0 +1,4 @@ +rhoPorousMRFPimpleFoam.C + +EXE = $(FOAM_APPBIN)/rhoPorousMRFPimpleFoam + diff --git a/applications/solvers/compressible/rhoPorousMRFPimpleFoam/Make/options b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/Make/options new file mode 100644 index 0000000000..61c1b6fe46 --- /dev/null +++ b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/Make/options @@ -0,0 +1,15 @@ +EXE_INC = \ + -I../rhoPimpleFoam \ + -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \ + -I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \ + -I$(LIB_SRC)/finiteVolume/cfdTools \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + -I$(LIB_SRC)/meshTools/lnInclude + +EXE_LIBS = \ + -lbasicThermophysicalModels \ + -lspecie \ + -lcompressibleRASModels \ + -lcompressibleLESModels \ + -lfiniteVolume \ + -lmeshTools diff --git a/applications/solvers/compressible/rhoPorousMRFPimpleFoam/UEqn.H b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/UEqn.H new file mode 100644 index 0000000000..3b4ead17e7 --- /dev/null +++ b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/UEqn.H @@ -0,0 +1,39 @@ +// Solve the Momentum equation + +tmp UEqn +( + pZones.ddt(rho, U) + + fvm::div(phi, U) + + turbulence->divDevRhoReff(U) +); + +if (oCorr == nOuterCorr-1) +{ + UEqn().relax(1); +} +else +{ + UEqn().relax(); +} + +mrfZones.addCoriolis(rho, UEqn()); +pZones.addResistance(UEqn()); + +volScalarField rUA = 1.0/UEqn().A(); + +if (momentumPredictor) +{ + if (oCorr == nOuterCorr-1) + { + solve(UEqn() == -fvc::grad(p), mesh.solver("UFinal")); + } + else + { + solve(UEqn() == -fvc::grad(p)); + } +} +else +{ + U = rUA*(UEqn().H() - fvc::grad(p)); + U.correctBoundaryConditions(); +} diff --git a/applications/solvers/compressible/rhoPorousMRFPimpleFoam/createFields.H b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/createFields.H new file mode 100644 index 0000000000..b9a86ef995 --- /dev/null +++ b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/createFields.H @@ -0,0 +1,70 @@ + Info<< "Reading thermophysical properties\n" << endl; + + autoPtr pThermo + ( + basicPsiThermo::New(mesh) + ); + basicPsiThermo& thermo = pThermo(); + + volScalarField& p = thermo.p(); + volScalarField& h = thermo.h(); + const volScalarField& psi = thermo.psi(); + + volScalarField rho + ( + IOobject + ( + "rho", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + thermo.rho() + ); + + Info<< "Reading field U\n" << endl; + volVectorField U + ( + IOobject + ( + "U", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + #include "compressibleCreatePhi.H" + + dimensionedScalar pMin + ( + mesh.solutionDict().subDict("PIMPLE").lookup("pMin") + ); + + Info<< "Creating turbulence model\n" << endl; + autoPtr turbulence + ( + compressible::turbulenceModel::New + ( + rho, + U, + phi, + thermo + ) + ); + + //dimensionedScalar initialMass = fvc::domainIntegrate(rho); + + + Info<< "Creating field DpDt\n" << endl; + volScalarField DpDt = + fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p); + + MRFZones mrfZones(mesh); + mrfZones.correctBoundaryVelocity(U); + + porousZones pZones(mesh); + Switch pressureImplicitPorosity(false); diff --git a/applications/solvers/compressible/rhoPorousMRFPimpleFoam/pEqn.H b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/pEqn.H new file mode 100644 index 0000000000..c74fb4d84b --- /dev/null +++ b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/pEqn.H @@ -0,0 +1,123 @@ +rho = thermo.rho(); + +volScalarField rUA = 1.0/UEqn().A(); +U = rUA*UEqn().H(); + +if (nCorr <= 1) +{ + UEqn.clear(); +} + +if (transonic) +{ + surfaceScalarField phid + ( + "phid", + fvc::interpolate(psi) + *( + (fvc::interpolate(U) & mesh.Sf()) + + fvc::ddtPhiCorr(rUA, rho, U, phi) + ) + ); + mrfZones.relativeFlux(fvc::interpolate(psi), phid); + + for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + { + fvScalarMatrix pEqn + ( + fvm::ddt(psi, p) + + fvm::div(phid, p) + - fvm::laplacian(rho*rUA, p) + ); + + if + ( + oCorr == nOuterCorr-1 + && corr == nCorr-1 + && nonOrth == nNonOrthCorr + ) + { + pEqn.solve(mesh.solver("pFinal")); + } + else + { + pEqn.solve(); + } + + if (nonOrth == nNonOrthCorr) + { + phi == pEqn.flux(); + } + } +} +else +{ + phi = + fvc::interpolate(rho)* + ( + (fvc::interpolate(U) & mesh.Sf()) + //+ fvc::ddtPhiCorr(rUA, rho, U, phi) + ); + mrfZones.relativeFlux(fvc::interpolate(rho), phi); + + for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + { + // Pressure corrector + fvScalarMatrix pEqn + ( + fvm::ddt(psi, p) + + fvc::div(phi) + - fvm::laplacian(rho*rUA, p) + ); + + if + ( + oCorr == nOuterCorr-1 + && corr == nCorr-1 + && nonOrth == nNonOrthCorr + ) + { + pEqn.solve(mesh.solver("pFinal")); + } + else + { + pEqn.solve(); + } + + if (nonOrth == nNonOrthCorr) + { + phi += pEqn.flux(); + } + } +} + +#include "rhoEqn.H" +#include "compressibleContinuityErrs.H" + +//if (oCorr != nOuterCorr-1) +{ + // Explicitly relax pressure for momentum corrector + p.relax(); + + rho = thermo.rho(); + rho.relax(); + Info<< "rho max/min : " << max(rho).value() + << " " << min(rho).value() << endl; +} + +U -= rUA*fvc::grad(p); +U.correctBoundaryConditions(); + +DpDt = fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p); + +bound(p, pMin); + +// For closed-volume cases adjust the pressure and density levels +// to obey overall mass continuity +/* +if (closedVolume) +{ + p += (initialMass - fvc::domainIntegrate(psi*p)) + /fvc::domainIntegrate(psi); +} +*/ diff --git a/applications/solvers/compressible/rhoPorousMRFPimpleFoam/rhoPorousMRFPimpleFoam.C b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/rhoPorousMRFPimpleFoam.C new file mode 100644 index 0000000000..45d7646d04 --- /dev/null +++ b/applications/solvers/compressible/rhoPorousMRFPimpleFoam/rhoPorousMRFPimpleFoam.C @@ -0,0 +1,105 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Application + rhoPorousMRFPimpleFoam + +Description + Transient solver for laminar or turbulent flow of compressible fluids + with support for porous media and MRF for HVAC and similar applications. + + Uses the flexible PIMPLE (PISO-SIMPLE) solution for time-resolved and + pseudo-transient simulations. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "basicPsiThermo.H" +#include "turbulenceModel.H" +#include "bound.H" +#include "MRFZones.H" +#include "porousZones.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + #include "setRootCase.H" + #include "createTime.H" + #include "createMesh.H" + #include "createFields.H" + #include "initContinuityErrs.H" + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + while (runTime.run()) + { + #include "readTimeControls.H" + #include "readPIMPLEControls.H" + #include "compressibleCourantNo.H" + #include "setDeltaT.H" + + runTime++; + + Info<< "Time = " << runTime.timeName() << nl << endl; + + if (nOuterCorr != 1) + { + p.storePrevIter(); + rho.storePrevIter(); + } + + #include "rhoEqn.H" + + // --- Pressure-velocity PIMPLE corrector loop + for (int oCorr=0; oCorrcorrect(); + } + + runTime.write(); + + Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s" + << " ClockTime = " << runTime.elapsedClockTime() << " s" + << nl << endl; + } + + Info<< "End\n" << endl; + + return 0; +} + + +// ************************************************************************* // diff --git a/applications/solvers/electromagnetics/mhdFoam/createFields.H b/applications/solvers/electromagnetics/mhdFoam/createFields.H index 5b346858b6..04e05c032d 100644 --- a/applications/solvers/electromagnetics/mhdFoam/createFields.H +++ b/applications/solvers/electromagnetics/mhdFoam/createFields.H @@ -61,7 +61,7 @@ mesh ); -# include "createPhi.H" + #include "createPhi.H" Info<< "Reading field pB\n" << endl; volScalarField pB @@ -93,15 +93,15 @@ ); -# include "createPhiB.H" + #include "createPhiB.H" -dimensionedScalar DB = 1.0/(mu*sigma); -DB.name() = "DB"; + dimensionedScalar DB = 1.0/(mu*sigma); + DB.name() = "DB"; -dimensionedScalar DBU = 1.0/(2.0*mu*rho); -DBU.name() = "DBU"; + dimensionedScalar DBU = 1.0/(2.0*mu*rho); + DBU.name() = "DBU"; -label pRefCell = 0; -scalar pRefValue = 0.0; -setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue); + label pRefCell = 0; + scalar pRefValue = 0.0; + setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue); diff --git a/applications/solvers/electromagnetics/mhdFoam/mhdFoam.C b/applications/solvers/electromagnetics/mhdFoam/mhdFoam.C index f09b35efc9..3603d04070 100644 --- a/applications/solvers/electromagnetics/mhdFoam/mhdFoam.C +++ b/applications/solvers/electromagnetics/mhdFoam/mhdFoam.C @@ -49,8 +49,6 @@ Description \*---------------------------------------------------------------------------*/ -#include "string.H" -#include "Time.H" #include "fvCFD.H" #include "OSspecific.H" diff --git a/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.C b/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.C index b1f8c3bed2..5013a4e6a2 100644 --- a/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.C +++ b/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.C @@ -38,7 +38,7 @@ solidWallHeatFluxTemperatureFvPatchScalarField const DimensionedField& iF ) : - fixedValueFvPatchScalarField(p, iF), + fixedGradientFvPatchScalarField(p, iF), q_(p.size(), 0.0), KName_("undefined-K") {} @@ -53,7 +53,7 @@ solidWallHeatFluxTemperatureFvPatchScalarField const fvPatchFieldMapper& mapper ) : - fixedValueFvPatchScalarField(ptf, p, iF, mapper), + fixedGradientFvPatchScalarField(ptf, p, iF, mapper), q_(ptf.q_, mapper), KName_(ptf.KName_) {} @@ -67,7 +67,7 @@ solidWallHeatFluxTemperatureFvPatchScalarField const dictionary& dict ) : - fixedValueFvPatchScalarField(p, iF, dict), + fixedGradientFvPatchScalarField(p, iF, dict), q_("q", dict, p.size()), KName_(dict.lookup("K")) {} @@ -79,7 +79,7 @@ solidWallHeatFluxTemperatureFvPatchScalarField const solidWallHeatFluxTemperatureFvPatchScalarField& tppsf ) : - fixedValueFvPatchScalarField(tppsf), + fixedGradientFvPatchScalarField(tppsf), q_(tppsf.q_), KName_(tppsf.KName_) {} @@ -92,7 +92,7 @@ solidWallHeatFluxTemperatureFvPatchScalarField const DimensionedField& iF ) : - fixedValueFvPatchScalarField(tppsf, iF), + fixedGradientFvPatchScalarField(tppsf, iF), q_(tppsf.q_), KName_(tppsf.KName_) {} @@ -105,7 +105,7 @@ void Foam::solidWallHeatFluxTemperatureFvPatchScalarField::autoMap const fvPatchFieldMapper& m ) { - fixedValueFvPatchScalarField::autoMap(m); + fixedGradientFvPatchScalarField::autoMap(m); q_.autoMap(m); } @@ -116,7 +116,7 @@ void Foam::solidWallHeatFluxTemperatureFvPatchScalarField::rmap const labelList& addr ) { - fixedValueFvPatchScalarField::rmap(ptf, addr); + fixedGradientFvPatchScalarField::rmap(ptf, addr); const solidWallHeatFluxTemperatureFvPatchScalarField& hfptf = refCast(ptf); @@ -132,14 +132,14 @@ void Foam::solidWallHeatFluxTemperatureFvPatchScalarField::updateCoeffs() return; } - const scalarField& Kw = - patch().lookupPatchField(KName_); + const scalarField& Kw = patch().lookupPatchField + ( + KName_ + ); - const fvPatchScalarField& Tw = *this; + gradient() = q_/Kw; - operator==(q_/(patch().deltaCoeffs()*Kw) + Tw.patchInternalField()); - - fixedValueFvPatchScalarField::updateCoeffs(); + fixedGradientFvPatchScalarField::updateCoeffs(); } @@ -148,9 +148,10 @@ void Foam::solidWallHeatFluxTemperatureFvPatchScalarField::write Ostream& os ) const { - fixedValueFvPatchScalarField::write(os); + fixedGradientFvPatchScalarField::write(os); q_.writeEntry("q", os); os.writeKeyword("K") << KName_ << token::END_STATEMENT << nl; + this->writeEntry("value", os); } diff --git a/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.H b/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.H index 85a1ef1cf1..d621e4d15a 100644 --- a/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.H +++ b/applications/solvers/heatTransfer/chtMultiRegionFoam/derivedFvPatchFields/solidWallHeatFluxTemperature/solidWallHeatFluxTemperatureFvPatchScalarField.H @@ -46,7 +46,7 @@ SourceFiles #ifndef solidWallHeatFluxTemperatureFvPatchScalarField_H #define solidWallHeatFluxTemperatureFvPatchScalarField_H -#include "fixedValueFvPatchFields.H" +#include "fixedGradientFvPatchFields.H" // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // @@ -59,7 +59,7 @@ namespace Foam class solidWallHeatFluxTemperatureFvPatchScalarField : - public fixedValueFvPatchScalarField + public fixedGradientFvPatchScalarField { // Private data diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/Make/files b/applications/solvers/incompressible/ajointShapeOptimizationFoam/Make/files new file mode 100644 index 0000000000..d42156d410 --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/Make/files @@ -0,0 +1,5 @@ +adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C +adjointOutletVelocity/adjointOutletVelocityFvPatchVectorField.C +adjointShapeOptimizationFoam.C + +EXE = $(FOAM_APPBIN)/adjointShapeOptimizationFoam diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/Make/options b/applications/solvers/incompressible/ajointShapeOptimizationFoam/Make/options new file mode 100644 index 0000000000..1223bdd06f --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/Make/options @@ -0,0 +1,11 @@ +EXE_INC = \ + -I$(LIB_SRC)/turbulenceModels \ + -I$(LIB_SRC)/turbulenceModels/incompressible/RAS/RASModel \ + -I$(LIB_SRC)/transportModels \ + -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \ + -I$(LIB_SRC)/finiteVolume/lnInclude + +EXE_LIBS = \ + -lincompressibleRASModels \ + -lincompressibleTransportModels \ + -lfiniteVolume diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointContinuityErrs.H b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointContinuityErrs.H new file mode 100644 index 0000000000..74ebd457c0 --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointContinuityErrs.H @@ -0,0 +1,47 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Global + continuityErrs + +Description + Calculates and prints the continuity errors. + +\*---------------------------------------------------------------------------*/ + +{ + scalar sumLocalContErr = runTime.deltaT().value()* + mag(fvc::div(phia))().weightedAverage(mesh.V()).value(); + + scalar globalContErr = runTime.deltaT().value()* + fvc::div(phia)().weightedAverage(mesh.V()).value(); + cumulativeContErr += globalContErr; + + Info<< "Adjoint continuity errors : sum local = " << sumLocalContErr + << ", global = " << globalContErr + << ", cumulative = " << cumulativeContErr + << endl; +} + +// ************************************************************************* // diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C new file mode 100644 index 0000000000..13debc1786 --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletPressure/adjointOutletPressureFvPatchScalarField.C @@ -0,0 +1,132 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "adjointOutletPressureFvPatchScalarField.H" +#include "addToRunTimeSelectionTable.H" +#include "fvPatchMapper.H" +#include "volFields.H" +#include "surfaceFields.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::adjointOutletPressureFvPatchScalarField:: +adjointOutletPressureFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF +) +: + fixedValueFvPatchScalarField(p, iF) +{} + + +Foam::adjointOutletPressureFvPatchScalarField:: +adjointOutletPressureFvPatchScalarField +( + const adjointOutletPressureFvPatchScalarField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + fixedValueFvPatchScalarField(ptf, p, iF, mapper) +{} + + +Foam::adjointOutletPressureFvPatchScalarField:: +adjointOutletPressureFvPatchScalarField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + fixedValueFvPatchScalarField(p, iF) +{ + fvPatchField::operator= + ( + scalarField("value", dict, p.size()) + ); +} + + +Foam::adjointOutletPressureFvPatchScalarField:: +adjointOutletPressureFvPatchScalarField +( + const adjointOutletPressureFvPatchScalarField& tppsf, + const DimensionedField& iF +) +: + fixedValueFvPatchScalarField(tppsf, iF) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +void Foam::adjointOutletPressureFvPatchScalarField::updateCoeffs() +{ + if (updated()) + { + return; + } + + const fvsPatchField& phip = + patch().lookupPatchField("phi"); + + const fvsPatchField& phiap = + patch().lookupPatchField("phia"); + + const fvPatchField& Up = + patch().lookupPatchField("U"); + + const fvPatchField& Uap = + patch().lookupPatchField("Ua"); + + operator==((phiap/patch().magSf() - 1.0)*phip/patch().magSf() + (Up & Uap)); + + fixedValueFvPatchScalarField::updateCoeffs(); +} + + +void Foam::adjointOutletPressureFvPatchScalarField::write(Ostream& os) const +{ + fvPatchScalarField::write(os); + writeEntry("value", os); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + makePatchTypeField + ( + fvPatchScalarField, + adjointOutletPressureFvPatchScalarField + ); +} + +// ************************************************************************* // diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletPressure/adjointOutletPressureFvPatchScalarField.H b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletPressure/adjointOutletPressureFvPatchScalarField.H new file mode 100644 index 0000000000..59514cb67f --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletPressure/adjointOutletPressureFvPatchScalarField.H @@ -0,0 +1,137 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Class + adjointOutletPressureFvPatchScalarField + +Description + +SourceFiles + adjointOutletPressureFvPatchScalarField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef adjointOutletPressureFvPatchScalarField_H +#define adjointOutletPressureFvPatchScalarField_H + +#include "fixedValueFvPatchFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class adjointOutletPressureFvPatch Declaration +\*---------------------------------------------------------------------------*/ + +class adjointOutletPressureFvPatchScalarField +: + public fixedValueFvPatchScalarField +{ + +public: + + //- Runtime type information + TypeName("adjointOutletPressure"); + + + // Constructors + + //- Construct from patch and internal field + adjointOutletPressureFvPatchScalarField + ( + const fvPatch&, + const DimensionedField& + ); + + //- Construct from patch, internal field and dictionary + adjointOutletPressureFvPatchScalarField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given adjointOutletPressureFvPatchScalarField + // onto a new patch + adjointOutletPressureFvPatchScalarField + ( + const adjointOutletPressureFvPatchScalarField&, + const fvPatch&, + const DimensionedField&, + const fvPatchFieldMapper& + ); + + //- Construct and return a clone + virtual tmp clone() const + { + return tmp + ( + new adjointOutletPressureFvPatchScalarField(*this) + ); + } + + //- Construct as copy setting internal field reference + adjointOutletPressureFvPatchScalarField + ( + const adjointOutletPressureFvPatchScalarField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp clone + ( + const DimensionedField& iF + ) const + { + return tmp + ( + new adjointOutletPressureFvPatchScalarField(*this, iF) + ); + } + + + // Member functions + + // Evaluation functions + + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + + //- Write + virtual void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletVelocity/adjointOutletVelocityFvPatchVectorField.C b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletVelocity/adjointOutletVelocityFvPatchVectorField.C new file mode 100644 index 0000000000..e99550c619 --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletVelocity/adjointOutletVelocityFvPatchVectorField.C @@ -0,0 +1,131 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "adjointOutletVelocityFvPatchVectorField.H" +#include "volFields.H" +#include "addToRunTimeSelectionTable.H" +#include "surfaceFields.H" +#include "fvPatchFieldMapper.H" + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::adjointOutletVelocityFvPatchVectorField:: +adjointOutletVelocityFvPatchVectorField +( + const fvPatch& p, + const DimensionedField& iF +) +: + fixedValueFvPatchVectorField(p, iF) +{} + + +Foam::adjointOutletVelocityFvPatchVectorField:: +adjointOutletVelocityFvPatchVectorField +( + const fvPatch& p, + const DimensionedField& iF, + const dictionary& dict +) +: + fixedValueFvPatchVectorField(p, iF) +{ + fvPatchVectorField::operator=(vectorField("value", dict, p.size())); +} + + +Foam::adjointOutletVelocityFvPatchVectorField:: +adjointOutletVelocityFvPatchVectorField +( + const adjointOutletVelocityFvPatchVectorField& ptf, + const fvPatch& p, + const DimensionedField& iF, + const fvPatchFieldMapper& mapper +) +: + fixedValueFvPatchVectorField(ptf, p, iF, mapper) +{} + + +Foam::adjointOutletVelocityFvPatchVectorField:: +adjointOutletVelocityFvPatchVectorField +( + const adjointOutletVelocityFvPatchVectorField& pivpvf, + const DimensionedField& iF +) +: + fixedValueFvPatchVectorField(pivpvf, iF) +{} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + +// Update the coefficients associated with the patch field +void Foam::adjointOutletVelocityFvPatchVectorField::updateCoeffs() +{ + if (updated()) + { + return; + } + + const fvsPatchField& phiap = + patch().lookupPatchField("phia"); + + const fvPatchField& Up = + patch().lookupPatchField("U"); + + scalarField Un = mag(patch().nf() & Up); + vectorField UtHat = (Up - patch().nf()*Un)/(Un + SMALL); + + vectorField Uan = patch().nf()*(patch().nf() & patchInternalField()); + + vectorField::operator=(phiap*patch().Sf()/sqr(patch().magSf()) + UtHat); + //vectorField::operator=(Uan + UtHat); + + fixedValueFvPatchVectorField::updateCoeffs(); +} + + +void Foam::adjointOutletVelocityFvPatchVectorField::write(Ostream& os) const +{ + fvPatchVectorField::write(os); + writeEntry("value", os); +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + makePatchTypeField + ( + fvPatchVectorField, + adjointOutletVelocityFvPatchVectorField + ); +} + + +// ************************************************************************* // diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletVelocity/adjointOutletVelocityFvPatchVectorField.H b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletVelocity/adjointOutletVelocityFvPatchVectorField.H new file mode 100644 index 0000000000..2e14f96f9b --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointOutletVelocity/adjointOutletVelocityFvPatchVectorField.H @@ -0,0 +1,131 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2007 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Class + adjointOutletVelocityFvPatchVectorField + +Description + +SourceFiles + adjointOutletVelocityFvPatchVectorField.C + +\*---------------------------------------------------------------------------*/ + +#ifndef adjointOutletVelocityFvPatchVectorField_H +#define adjointOutletVelocityFvPatchVectorField_H + +#include "fixedValueFvPatchFields.H" + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +namespace Foam +{ + +/*---------------------------------------------------------------------------*\ + Class adjointOutletVelocityFvPatch Declaration +\*---------------------------------------------------------------------------*/ + +class adjointOutletVelocityFvPatchVectorField +: + public fixedValueFvPatchVectorField +{ + +public: + + //- Runtime type information + TypeName("adjointOutletVelocity"); + + + // Constructors + + //- Construct from patch and internal field + adjointOutletVelocityFvPatchVectorField + ( + const fvPatch&, + const DimensionedField& + ); + + //- Construct from patch, internal field and dictionary + adjointOutletVelocityFvPatchVectorField + ( + const fvPatch&, + const DimensionedField&, + const dictionary& + ); + + //- Construct by mapping given adjointOutletVelocityFvPatchVectorField + // onto a new patch + adjointOutletVelocityFvPatchVectorField + ( + const adjointOutletVelocityFvPatchVectorField&, + const fvPatch&, + const DimensionedField&, + const fvPatchFieldMapper& + ); + + //- Construct and return a clone + virtual tmp clone() const + { + return tmp + ( + new adjointOutletVelocityFvPatchVectorField(*this) + ); + } + + //- Construct as copy setting internal field reference + adjointOutletVelocityFvPatchVectorField + ( + const adjointOutletVelocityFvPatchVectorField&, + const DimensionedField& + ); + + //- Construct and return a clone setting internal field reference + virtual tmp clone(const DimensionedField& iF) const + { + return tmp + ( + new adjointOutletVelocityFvPatchVectorField(*this, iF) + ); + } + + + // Member functions + + //- Update the coefficients associated with the patch field + virtual void updateCoeffs(); + + //- Write + virtual void write(Ostream&) const; +}; + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +} // End namespace Foam + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointShapeOptimizationFoam.C b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointShapeOptimizationFoam.C new file mode 100644 index 0000000000..acf1be057f --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/adjointShapeOptimizationFoam.C @@ -0,0 +1,226 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Application + ajointShapeOptimizationFoam + +Description + Steady-state solver for incompressible, turbulent flow of non-Newtonian + fluids with optimisation of duct shape by applying "blockage" in regions + causing pressure loss as estimated using an adjoint formulation. + + References: + @verbatim + "Implementation of a continuous adjoint for topology optimization of + ducted flows" + C. Othmer, + E. de Villiers, + H.G. Weller + AIAA-2007-3947 + http://pdf.aiaa.org/preview/CDReadyMCFD07_1379/PV2007_3947.pdf + @endverbatim + + Note that this solver optimises for total pressure loss whereas the + above paper describes the method for optimising power-loss. + +\*---------------------------------------------------------------------------*/ + +#include "fvCFD.H" +#include "singlePhaseTransportModel.H" +#include "RASModel.H" + +template +void zeroCells +( + GeometricField& vf, + const labelList& cells +) +{ + forAll(cells, i) + { + vf[cells[i]] = pTraits::zero; + } +} + + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +int main(int argc, char *argv[]) +{ + #include "setRootCase.H" + + #include "createTime.H" + #include "createMesh.H" + #include "createFields.H" + #include "initContinuityErrs.H" + + // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + + Info<< "\nStarting time loop\n" << endl; + + for (runTime++; !runTime.end(); runTime++) + { + Info<< "Time = " << runTime.timeName() << nl << endl; + + #include "readSIMPLEControls.H" + + p.storePrevIter(); + + laminarTransport.lookup("lambda") >> lambda; + + //alpha += + // mesh.relaxationFactor("alpha") + // *(lambda*max(Ua & U, zeroSensitivity) - alpha); + alpha += + mesh.relaxationFactor("alpha") + *(min(max(alpha + lambda*(Ua & U), zeroAlpha), alphaMax) - alpha); + + zeroCells(alpha, inletCells); + //zeroCells(alpha, outletCells); + + // Pressure-velocity SIMPLE corrector + { + // Momentum predictor + + tmp UEqn + ( + fvm::div(phi, U) + + turbulence->divDevReff(U) + + fvm::Sp(alpha, U) + ); + + UEqn().relax(); + + solve(UEqn() == -fvc::grad(p)); + + p.boundaryField().updateCoeffs(); + volScalarField rAU = 1.0/UEqn().A(); + U = rAU*UEqn().H(); + UEqn.clear(); + phi = fvc::interpolate(U) & mesh.Sf(); + adjustPhi(phi, U, p); + + // Non-orthogonal pressure corrector loop + for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + { + fvScalarMatrix pEqn + ( + fvm::laplacian(rAU, p) == fvc::div(phi) + ); + + pEqn.setReference(pRefCell, pRefValue); + pEqn.solve(); + + if (nonOrth == nNonOrthCorr) + { + phi -= pEqn.flux(); + } + } + + #include "continuityErrs.H" + + // Explicitly relax pressure for momentum corrector + p.relax(); + + // Momentum corrector + U -= rAU*fvc::grad(p); + U.correctBoundaryConditions(); + } + + pa.storePrevIter(); + + // Adjoint Pressure-velocity SIMPLE corrector + { + // Adjoint Momentum predictor + + volVectorField adjointTransposeConvection = (fvc::grad(Ua) & U); + //volVectorField adjointTransposeConvection = fvc::reconstruct + //( + // mesh.magSf()*(fvc::snGrad(Ua) & fvc::interpolate(U)) + //); + + zeroCells(adjointTransposeConvection, inletCells); + + tmp UaEqn + ( + fvm::div(-phi, Ua) + - adjointTransposeConvection + + turbulence->divDevReff(Ua) + + fvm::Sp(alpha, Ua) + ); + + UaEqn().relax(); + + solve(UaEqn() == -fvc::grad(pa)); + + pa.boundaryField().updateCoeffs(); + volScalarField rAUa = 1.0/UaEqn().A(); + Ua = rAUa*UaEqn().H(); + UaEqn.clear(); + phia = fvc::interpolate(Ua) & mesh.Sf(); + adjustPhi(phia, Ua, pa); + + // Non-orthogonal pressure corrector loop + for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) + { + fvScalarMatrix paEqn + ( + fvm::laplacian(rAUa, pa) == fvc::div(phia) + ); + + paEqn.setReference(paRefCell, paRefValue); + paEqn.solve(); + + if (nonOrth == nNonOrthCorr) + { + phia -= paEqn.flux(); + } + } + + #include "adjointContinuityErrs.H" + + // Explicitly relax pressure for adjoint momentum corrector + pa.relax(); + + // Adjoint momentum corrector + Ua -= rAUa*fvc::grad(pa); + Ua.correctBoundaryConditions(); + } + + turbulence->correct(); + + runTime.write(); + + Info<< "ExecutionTime = " + << runTime.elapsedCpuTime() + << " s\n\n" << endl; + } + + Info<< "End\n" << endl; + + return(0); +} + + +// ************************************************************************* // diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/createFields.H b/applications/solvers/incompressible/ajointShapeOptimizationFoam/createFields.H new file mode 100644 index 0000000000..dd2f6416ee --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/createFields.H @@ -0,0 +1,105 @@ + Info << "Reading field p\n" << endl; + volScalarField p + ( + IOobject + ( + "p", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + Info << "Reading field U\n" << endl; + volVectorField U + ( + IOobject + ( + "U", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + #include "createPhi.H" + + + label pRefCell = 0; + scalar pRefValue = 0.0; + setRefCell(p, mesh.solutionDict().subDict("SIMPLE"), pRefCell, pRefValue); + + + Info << "Reading field pa\n" << endl; + volScalarField pa + ( + IOobject + ( + "pa", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + Info << "Reading field Ua\n" << endl; + volVectorField Ua + ( + IOobject + ( + "Ua", + runTime.timeName(), + mesh, + IOobject::MUST_READ, + IOobject::AUTO_WRITE + ), + mesh + ); + + #include "createPhia.H" + + + label paRefCell = 0; + scalar paRefValue = 0.0; + setRefCell(pa, mesh.solutionDict().subDict("SIMPLE"), paRefCell, paRefValue); + + + singlePhaseTransportModel laminarTransport(U, phi); + + autoPtr turbulence + ( + incompressible::RASModel::New(U, phi, laminarTransport) + ); + + + dimensionedScalar zeroSensitivity("0", dimVelocity*dimVelocity, 0.0); + dimensionedScalar zeroAlpha("0", dimless/dimTime, 0.0); + + dimensionedScalar lambda(laminarTransport.lookup("lambda")); + dimensionedScalar alphaMax(laminarTransport.lookup("alphaMax")); + + const labelList& inletCells = + mesh.boundary()[mesh.boundaryMesh().findPatchID("inlet")].faceCells(); + //const labelList& outletCells = + // mesh.boundary()[mesh.boundaryMesh().findPatchID("outlet")].faceCells(); + + volScalarField alpha + ( + IOobject + ( + "alpha", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + lambda*max(Ua & U, zeroSensitivity) + ); + zeroCells(alpha, inletCells); + //zeroCells(alpha, outletCells); diff --git a/applications/solvers/incompressible/ajointShapeOptimizationFoam/createPhia.H b/applications/solvers/incompressible/ajointShapeOptimizationFoam/createPhia.H new file mode 100644 index 0000000000..ea227ab133 --- /dev/null +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/createPhia.H @@ -0,0 +1,57 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +Global + createPhia + +Description + Creates and initialises the face-flux field phia. + +\*---------------------------------------------------------------------------*/ + +#ifndef createPhia_H +#define createPhia_H + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +Info<< "Reading/calculating face flux field phia\n" << endl; + +surfaceScalarField phia +( + IOobject + ( + "phia", + runTime.timeName(), + mesh, + IOobject::READ_IF_PRESENT, + IOobject::AUTO_WRITE + ), + linearInterpolate(Ua) & mesh.Sf() +); + +// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // + +#endif + +// ************************************************************************* // diff --git a/applications/solvers/multiphase/compressibleInterDyMFoam/alphaEqns.H b/applications/solvers/multiphase/compressibleInterDyMFoam/alphaEqns.H index 819cd0f538..f545262767 100644 --- a/applications/solvers/multiphase/compressibleInterDyMFoam/alphaEqns.H +++ b/applications/solvers/multiphase/compressibleInterDyMFoam/alphaEqns.H @@ -59,7 +59,17 @@ alpharScheme ); - MULES::explicitSolve(oneField(), alpha1, phi, phiAlpha1, Sp, Su, 1, 0); + MULES::explicitSolve + ( + geometricOneField(), + alpha1, + phi, + phiAlpha1, + Sp, + Su, + 1, + 0 + ); surfaceScalarField rho1f = fvc::interpolate(rho1); surfaceScalarField rho2f = fvc::interpolate(rho2); diff --git a/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H b/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H index 819cd0f538..dd704c0693 100644 --- a/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H +++ b/applications/solvers/multiphase/compressibleInterFoam/alphaEqns.H @@ -59,7 +59,7 @@ alpharScheme ); - MULES::explicitSolve(oneField(), alpha1, phi, phiAlpha1, Sp, Su, 1, 0); + MULES::explicitSolve(geometricOneField(), alpha1, phi, phiAlpha1, Sp, Su, 1, 0); surfaceScalarField rho1f = fvc::interpolate(rho1); surfaceScalarField rho2f = fvc::interpolate(rho2); diff --git a/applications/solvers/multiphase/interMixingFoam/alphaEqns.H b/applications/solvers/multiphase/interMixingFoam/alphaEqns.H index f9bad0a705..972d1d7121 100644 --- a/applications/solvers/multiphase/interMixingFoam/alphaEqns.H +++ b/applications/solvers/multiphase/interMixingFoam/alphaEqns.H @@ -70,12 +70,12 @@ MULES::limiter ( allLambda, - oneField(), + geometricOneField(), alpha1, phiAlpha1BD, phiAlpha1, - zeroField(), - zeroField(), + zero(), + zero(), 1, 0, 3 @@ -107,12 +107,12 @@ MULES::limiter ( allLambda, - oneField(), + geometricOneField(), alpha2, phiAlpha2BD, phiAlpha2, - zeroField(), - zeroField(), + zero(), + zero(), 1, 0, 3 diff --git a/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H b/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H index 15a5291ee6..a232056b8d 100644 --- a/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H +++ b/applications/solvers/multiphase/interPhaseChangeFoam/alphaEqn.H @@ -51,8 +51,8 @@ ); //MULES::explicitSolve(alpha1, phi, phiAlpha, 1, 0); - //MULES::explicitSolve(oneField(), alpha1, phi, phiAlpha, Sp, Su, 1, 0); - MULES::implicitSolve(oneField(), alpha1, phi, phiAlpha, Sp, Su, 1, 0); + //MULES::explicitSolve(geometricOneField(), alpha1, phi, phiAlpha, Sp, Su, 1, 0); + MULES::implicitSolve(geometricOneField(), alpha1, phi, phiAlpha, Sp, Su, 1, 0); rhoPhi += (runTime.deltaT()/totalDeltaT) diff --git a/applications/test/nearWallDist-wave/testWallDistData.C b/applications/test/nearWallDist-wave/testWallDistData.C index cb742f435f..0b94dfa216 100644 --- a/applications/test/nearWallDist-wave/testWallDistData.C +++ b/applications/test/nearWallDist-wave/testWallDistData.C @@ -58,7 +58,7 @@ int main(int argc, char *argv[]) mesh ), mesh, - dimensionedVector("n", dimLength, wallPoint::greatPoint) + dimensionedVector("n", dimLength, point::max) ); // Fill wall patches with unit normal diff --git a/applications/test/syncTools/syncToolsTest.C b/applications/test/syncTools/syncToolsTest.C index 320b4e6f61..fb225ac397 100644 --- a/applications/test/syncTools/syncToolsTest.C +++ b/applications/test/syncTools/syncToolsTest.C @@ -187,8 +187,6 @@ void testSparseData(const polyMesh& mesh, Random& rndGen) ); const pointField& localPoints = allBoundary.localPoints(); - const point greatPoint(GREAT, GREAT, GREAT); - // Point data // ~~~~~~~~~~ @@ -196,7 +194,7 @@ void testSparseData(const polyMesh& mesh, Random& rndGen) { // Create some data. Use slightly perturbed positions. Map sparseData; - pointField fullData(mesh.nPoints(), greatPoint); + pointField fullData(mesh.nPoints(), point::max); forAll(localPoints, i) { @@ -222,7 +220,7 @@ void testSparseData(const polyMesh& mesh, Random& rndGen) mesh, fullData, minEqOp(), - greatPoint, + point::max, true // apply separation ); @@ -232,7 +230,7 @@ void testSparseData(const polyMesh& mesh, Random& rndGen) { const point& fullPt = fullData[meshPointI]; - if (fullPt != greatPoint) + if (fullPt != point::max) { const point& sparsePt = sparseData[meshPointI]; @@ -272,7 +270,7 @@ void testSparseData(const polyMesh& mesh, Random& rndGen) { // Create some data. Use slightly perturbed positions. EdgeMap sparseData; - pointField fullData(mesh.nEdges(), greatPoint); + pointField fullData(mesh.nEdges(), point::max); const edgeList& edges = allBoundary.edges(); const labelList meshEdges = allBoundary.meshEdges @@ -307,7 +305,7 @@ void testSparseData(const polyMesh& mesh, Random& rndGen) mesh, fullData, minEqOp(), - greatPoint, + point::max, true ); @@ -317,7 +315,7 @@ void testSparseData(const polyMesh& mesh, Random& rndGen) { const point& fullPt = fullData[meshEdgeI]; - if (fullPt != greatPoint) + if (fullPt != point::max) { const point& sparsePt = sparseData[mesh.edges()[meshEdgeI]]; @@ -364,8 +362,6 @@ void testPointSync(const polyMesh& mesh, Random& rndGen) { Info<< nl << "Testing point-wise data synchronisation." << endl; - const point greatPoint(GREAT, GREAT, GREAT); - // Test position. { @@ -379,7 +375,7 @@ void testPointSync(const polyMesh& mesh, Random& rndGen) mesh, syncedPoints, minEqOp(), - greatPoint, + point::max, true ); @@ -444,8 +440,6 @@ void testEdgeSync(const polyMesh& mesh, Random& rndGen) const edgeList& edges = mesh.edges(); - const point greatPoint(GREAT, GREAT, GREAT); - // Test position. { @@ -463,7 +457,7 @@ void testEdgeSync(const polyMesh& mesh, Random& rndGen) mesh, syncedMids, minEqOp(), - greatPoint, + point::max, true ); diff --git a/applications/utilities/mesh/generation/extrudeToRegionMesh/Make/files b/applications/utilities/mesh/generation/extrudeToRegionMesh/Make/files new file mode 100644 index 0000000000..62431e23d4 --- /dev/null +++ b/applications/utilities/mesh/generation/extrudeToRegionMesh/Make/files @@ -0,0 +1,7 @@ +patchPointEdgeCirculator.C +createShellMesh.C +extrudeToRegionMesh.C + +EXE = $(FOAM_APPBIN)/extrudeToRegionMesh + + diff --git a/applications/utilities/mesh/generation/extrudeToRegionMesh/Make/options b/applications/utilities/mesh/generation/extrudeToRegionMesh/Make/options new file mode 100644 index 0000000000..dbf50a9b09 --- /dev/null +++ b/applications/utilities/mesh/generation/extrudeToRegionMesh/Make/options @@ -0,0 +1,14 @@ +EXE_INC = \ + -I$(LIB_SRC)/finiteVolume/lnInclude \ + /* -I$(LIB_SRC)/surfMesh/lnInclude */ \ + /* -I$(LIB_SRC)/lagrangian/basic/lnInclude */ \ + -I$(LIB_SRC)/meshTools/lnInclude \ + -I$(LIB_SRC)/dynamicMesh/lnInclude + +EXE_LIBS = \ + -lfiniteVolume \ + /* -lsurfMesh */ \ + /* -llagrangian */ \ + -lmeshTools \ + -ldynamicMesh + diff --git a/applications/utilities/mesh/generation/extrudeToRegionMesh/createShellMesh.C b/applications/utilities/mesh/generation/extrudeToRegionMesh/createShellMesh.C new file mode 100644 index 0000000000..df998941ef --- /dev/null +++ b/applications/utilities/mesh/generation/extrudeToRegionMesh/createShellMesh.C @@ -0,0 +1,588 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / F ield | OpenFOAM: The Open Source CFD Toolbox + \\ / O peration | + \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd. + \\/ M anipulation | +------------------------------------------------------------------------------- +License + This file is part of OpenFOAM. + + OpenFOAM is free software; you can redistribute it and/or modify it + under the terms of the GNU General Public License as published by the + Free Software Foundation; either version 2 of the License, or (at your + option) any later version. + + OpenFOAM is distributed in the hope that it will be useful, but WITHOUT + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + for more details. + + You should have received a copy of the GNU General Public License + along with OpenFOAM; if not, write to the Free Software Foundation, + Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA + +\*---------------------------------------------------------------------------*/ + +#include "createShellMesh.H" +#include "polyTopoChange.H" +#include "meshTools.H" +#include "mapPolyMesh.H" +#include "polyAddPoint.H" +#include "polyAddFace.H" +#include "polyModifyFace.H" +#include "polyAddCell.H" +#include "patchPointEdgeCirculator.H" + +// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * // + +defineTypeNameAndDebug(Foam::createShellMesh, 0); + + +// * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * // + +void Foam::createShellMesh::calcPointRegions +( + const primitiveFacePatch& patch, + const PackedBoolList& nonManifoldEdge, + faceList& pointRegions, + labelList& regionPoints +) +{ + pointRegions.setSize(patch.size()); + forAll(pointRegions, faceI) + { + const face& f = patch.localFaces()[faceI]; + pointRegions[faceI].setSize(f.size(), -1); + } + + label nRegions = 0; + + forAll(pointRegions, faceI) + { + const face& f = patch.localFaces()[faceI]; + + forAll(f, fp) + { + if (pointRegions[faceI][fp] == -1) + { + // Found unassigned point. Distribute current region. + label pointI = f[fp]; + label edgeI = patch.faceEdges()[faceI][fp]; + + patchPointEdgeCirculator circ + ( + patch, + nonManifoldEdge, + edgeI, + findIndex(patch.edgeFaces()[edgeI], faceI), + pointI + ); + + for + ( + patchPointEdgeCirculator iter = circ.begin(); + iter != circ.end(); + ++iter + ) + { + label face2 = iter.faceID(); + + if (face2 != -1) + { + const face& f2 = patch.localFaces()[face2]; + label fp2 = findIndex(f2, pointI); + label& region = pointRegions[face2][fp2]; + if (region != -1) + { + FatalErrorIn + ( + "createShellMesh::calcPointRegions(..)" + ) << "On point " << pointI + << " at:" << patch.localPoints()[pointI] + << " found region:" << region + << abort(FatalError); + } + region = nRegions; + } + } + + nRegions++; + } + } + } + + + // From region back to originating point (many to one, a point might + // have multiple regions though) + regionPoints.setSize(nRegions); + forAll(pointRegions, faceI) + { + const face& f = patch.localFaces()[faceI]; + + forAll(f, fp) + { + regionPoints[pointRegions[faceI][fp]] = f[fp]; + } + } + + + if (debug) + { + const labelListList& pointFaces = patch.pointFaces(); + forAll(pointFaces, pointI) + { + label region = -1; + const labelList& pFaces = pointFaces[pointI]; + forAll(pFaces, i) + { + label faceI = pFaces[i]; + const face& f = patch.localFaces()[faceI]; + label fp = findIndex(f, pointI); + + if (region == -1) + { + region = pointRegions[faceI][fp]; + } + else if (region != pointRegions[faceI][fp]) + { + Pout<< "Non-manifold point:" << pointI + << " at " << patch.localPoints()[pointI] + << " region:" << region + << " otherRegion:" << pointRegions[faceI][fp] + << endl; + + } + } + } + } +} + + + +// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * // + +Foam::createShellMesh::createShellMesh +( + const primitiveFacePatch& patch, + const faceList& pointRegions, + const labelList& regionPoints +) +: + patch_(patch), + pointRegions_(pointRegions), + regionPoints_(regionPoints) +{ + if (pointRegions_.size() != patch_.size()) + { + FatalErrorIn("createShellMesh::createShellMesh(..)") + << "nFaces:" << patch_.size() + << " pointRegions:" << pointRegions.size() + << exit(FatalError); + } +} + + +// * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * // + + +void Foam::createShellMesh::setRefinement +( + const pointField& thickness, + const labelList& topPatchID, + const labelList& bottomPatchID, + const labelListList& extrudeEdgePatches, + polyTopoChange& meshMod +) +{ + if (thickness.size() != regionPoints_.size()) + { + FatalErrorIn("createShellMesh::setRefinement(..)") + << "nRegions:" << regionPoints_.size() + << " thickness:" << thickness.size() + << exit(FatalError); + } + + if + ( + topPatchID.size() != patch_.size() + && bottomPatchID.size() != patch_.size() + ) + { + FatalErrorIn("createShellMesh::setRefinement(..)") + << "nFaces:" << patch_.size() + << " topPatchID:" << topPatchID.size() + << " bottomPatchID:" << bottomPatchID.size() + << exit(FatalError); + } + + if (extrudeEdgePatches.size() != patch_.nEdges()) + { + FatalErrorIn("createShellMesh::setRefinement(..)") + << "nEdges:" << patch_.nEdges() + << " extrudeEdgePatches:" << extrudeEdgePatches.size() + << exit(FatalError); + } + + + + // From cell to patch (trivial) + DynamicList