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..e91cbf2ed6 100644 --- a/README +++ b/README @@ -112,7 +112,7 @@ 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/platform/$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. @@ -166,6 +166,7 @@ gcc-4.3.3. Execute the following: + cd $WM_THIRD_PARTY_DIR + rm -rf paraview-3.6.1/platforms + + rm -rf platforms/*/paraview-3.6.1 + ./makeParaView The PV3blockMeshReader and the PV3FoamReader ParaView plugins are compiled 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/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/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/MC/MC.C b/applications/solvers/incompressible/ajointShapeOptimizationFoam/createPhia.H similarity index 70% rename from src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/MC/MC.C rename to applications/solvers/incompressible/ajointShapeOptimizationFoam/createPhia.H index 901d524a8c..ea227ab133 100644 --- a/src/finiteVolume/interpolation/surfaceInterpolation/limitedSchemes/MC/MC.C +++ b/applications/solvers/incompressible/ajointShapeOptimizationFoam/createPhia.H @@ -2,7 +2,7 @@ ========= | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox \\ / O peration | - \\ / A nd | Copyright (C) 1991-2009 OpenCFD Ltd. + \\ / A nd | Copyright (C) 1991-2010 OpenCFD Ltd. \\/ M anipulation | ------------------------------------------------------------------------------- License @@ -22,17 +22,36 @@ 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. + \*---------------------------------------------------------------------------*/ -#include "LimitedScheme.H" -#include "MC.H" +#ifndef createPhia_H +#define createPhia_H // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // -namespace Foam -{ - makeLimitedSurfaceInterpolationScheme(MC, MCLimiter) - makeLimitedVSurfaceInterpolationScheme(MCV, MCLimiter) -} +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/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..da1fc43066 --- /dev/null +++ b/applications/utilities/mesh/generation/extrudeToRegionMesh/createShellMesh.C @@ -0,0 +1,592 @@ +/*---------------------------------------------------------------------------*\ + ========= | + \\ / 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& extrudeMasterPatchID, + const labelList& extrudeSlavePatchID, + const labelListList& extrudeEdgePatches, + polyTopoChange& meshMod +) +{ + if (thickness.size() != regionPoints_.size()) + { + FatalErrorIn("createShellMesh::setRefinement(..)") + << "nRegions:" << regionPoints_.size() + << " thickness:" << thickness.size() + << exit(FatalError); + } + + if + ( + extrudeMasterPatchID.size() != patch_.size() + && extrudeSlavePatchID.size() != patch_.size() + ) + { + FatalErrorIn("createShellMesh::setRefinement(..)") + << "nFaces:" << patch_.size() + << " extrudeMasterPatchID:" << extrudeMasterPatchID.size() + << " extrudeSlavePatchID:" << extrudeSlavePatchID.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