00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026 #include "tractionDisplacementFvPatchVectorField.H"
00027 #include <OpenFOAM/addToRunTimeSelectionTable.H>
00028 #include <finiteVolume/volFields.H>
00029
00030
00031
00032 namespace Foam
00033 {
00034
00035
00036
00037 tractionDisplacementFvPatchVectorField::
00038 tractionDisplacementFvPatchVectorField
00039 (
00040 const fvPatch& p,
00041 const DimensionedField<vector, volMesh>& iF
00042 )
00043 :
00044 fixedGradientFvPatchVectorField(p, iF),
00045 traction_(p.size(), vector::zero),
00046 pressure_(p.size(), 0.0)
00047 {
00048 fvPatchVectorField::operator=(patchInternalField());
00049 gradient() = vector::zero;
00050 }
00051
00052
00053 tractionDisplacementFvPatchVectorField::
00054 tractionDisplacementFvPatchVectorField
00055 (
00056 const tractionDisplacementFvPatchVectorField& tdpvf,
00057 const fvPatch& p,
00058 const DimensionedField<vector, volMesh>& iF,
00059 const fvPatchFieldMapper& mapper
00060 )
00061 :
00062 fixedGradientFvPatchVectorField(tdpvf, p, iF, mapper),
00063 traction_(tdpvf.traction_, mapper),
00064 pressure_(tdpvf.pressure_, mapper)
00065 {}
00066
00067
00068 tractionDisplacementFvPatchVectorField::
00069 tractionDisplacementFvPatchVectorField
00070 (
00071 const fvPatch& p,
00072 const DimensionedField<vector, volMesh>& iF,
00073 const dictionary& dict
00074 )
00075 :
00076 fixedGradientFvPatchVectorField(p, iF),
00077 traction_("traction", dict, p.size()),
00078 pressure_("pressure", dict, p.size())
00079 {
00080 fvPatchVectorField::operator=(patchInternalField());
00081 gradient() = vector::zero;
00082 }
00083
00084
00085 tractionDisplacementFvPatchVectorField::
00086 tractionDisplacementFvPatchVectorField
00087 (
00088 const tractionDisplacementFvPatchVectorField& tdpvf
00089 )
00090 :
00091 fixedGradientFvPatchVectorField(tdpvf),
00092 traction_(tdpvf.traction_),
00093 pressure_(tdpvf.pressure_)
00094 {}
00095
00096
00097 tractionDisplacementFvPatchVectorField::
00098 tractionDisplacementFvPatchVectorField
00099 (
00100 const tractionDisplacementFvPatchVectorField& tdpvf,
00101 const DimensionedField<vector, volMesh>& iF
00102 )
00103 :
00104 fixedGradientFvPatchVectorField(tdpvf, iF),
00105 traction_(tdpvf.traction_),
00106 pressure_(tdpvf.pressure_)
00107 {}
00108
00109
00110
00111
00112 void tractionDisplacementFvPatchVectorField::autoMap
00113 (
00114 const fvPatchFieldMapper& m
00115 )
00116 {
00117 fixedGradientFvPatchVectorField::autoMap(m);
00118 traction_.autoMap(m);
00119 pressure_.autoMap(m);
00120 }
00121
00122
00123 void tractionDisplacementFvPatchVectorField::rmap
00124 (
00125 const fvPatchVectorField& ptf,
00126 const labelList& addr
00127 )
00128 {
00129 fixedGradientFvPatchVectorField::rmap(ptf, addr);
00130
00131 const tractionDisplacementFvPatchVectorField& dmptf =
00132 refCast<const tractionDisplacementFvPatchVectorField>(ptf);
00133
00134 traction_.rmap(dmptf.traction_, addr);
00135 pressure_.rmap(dmptf.pressure_, addr);
00136 }
00137
00138
00139 void tractionDisplacementFvPatchVectorField::updateCoeffs()
00140 {
00141 if (updated())
00142 {
00143 return;
00144 }
00145
00146 const dictionary& mechanicalProperties =
00147 db().lookupObject<IOdictionary>("mechanicalProperties");
00148
00149 const dictionary& thermalProperties =
00150 db().lookupObject<IOdictionary>("thermalProperties");
00151
00152 dimensionedScalar rho(mechanicalProperties.lookup("rho"));
00153 dimensionedScalar rhoE(mechanicalProperties.lookup("E"));
00154 dimensionedScalar nu(mechanicalProperties.lookup("nu"));
00155
00156 dimensionedScalar E = rhoE/rho;
00157 dimensionedScalar mu = E/(2.0*(1.0 + nu));
00158 dimensionedScalar lambda = nu*E/((1.0 + nu)*(1.0 - 2.0*nu));
00159 dimensionedScalar threeK = E/(1.0 - 2.0*nu);
00160
00161 Switch planeStress(mechanicalProperties.lookup("planeStress"));
00162
00163 if (planeStress)
00164 {
00165 lambda = nu*E/((1.0 + nu)*(1.0 - nu));
00166 threeK = E/(1.0 - nu);
00167 }
00168
00169 scalar twoMuLambda = (2*mu + lambda).value();
00170
00171 vectorField n = patch().nf();
00172
00173 const fvPatchField<symmTensor>& sigmaD =
00174 patch().lookupPatchField<volSymmTensorField, symmTensor>("sigmaD");
00175
00176 gradient() =
00177 (
00178 (traction_ + pressure_*n)/rho.value()
00179 + twoMuLambda*fvPatchField<vector>::snGrad() - (n & sigmaD)
00180 )/twoMuLambda;
00181
00182 Switch thermalStress(thermalProperties.lookup("thermalStress"));
00183
00184 if (thermalStress)
00185 {
00186 dimensionedScalar alpha(thermalProperties.lookup("alpha"));
00187 dimensionedScalar threeKalpha = threeK*alpha;
00188
00189 const fvPatchField<scalar>& T =
00190 patch().lookupPatchField<volScalarField, scalar>("T");
00191
00192 gradient() += n*threeKalpha.value()*T/twoMuLambda;
00193 }
00194
00195 fixedGradientFvPatchVectorField::updateCoeffs();
00196 }
00197
00198
00199 void tractionDisplacementFvPatchVectorField::write(Ostream& os) const
00200 {
00201 fvPatchVectorField::write(os);
00202 traction_.writeEntry("traction", os);
00203 pressure_.writeEntry("pressure", os);
00204 writeEntry("value", os);
00205 }
00206
00207
00208
00209
00210 makePatchTypeField(fvPatchVectorField, tractionDisplacementFvPatchVectorField);
00211
00212
00213
00214 }
00215
00216