Commit e6dfcdf1 authored by Valentin Vallaeys's avatar Valentin Vallaeys
Browse files

rm boundary WallLinear and OpenFree for mom3D

parent 113a8aa0
......@@ -447,39 +447,6 @@ dgBoundaryCondition *dgConservationLawSW3dMomentum::newBoundaryWallVerticalBotto
return new dgConservationLawSW3dMomentumBoundaryWallVerticalBottom(this, extDiffFlux, zBotDist);
}
void dgConservationLawSW3dMomentum::boundaryWallLinear(functorCache &cm, fullMatrix<double> &val, double uvPenaltyFactor) const {
const fullMatrix<double> &normals = cm.get(*function::getNormals());
const fullMatrix<double> &uv = cm.get(*function::getSolution());
const fullMatrix<double> &uvGrad = cm.get(*function::getSolutionGradient());
const fullMatrix<double> *nuH = _nuH ? &cm.get(*_nuH) : NULL;
// bnd condition for linear momentum
val.resize(cm.nEvaluationPoint(), 2, true); //already set to zero in resize
for(int i=0;i<cm.nEvaluationPoint(); i++) {
double nx = normals(i,0);
double ny = normals(i,1);
double nz = normals(i,2);
double un = uv(i,0)*nx + uv(i,1)*ny;
if ( fabs(nz) < 1e-8 ) { // only for lateral faces!
val(i,0) = 0;
val(i,1) = 0;
val(i,0) += -uvPenaltyFactor*un*nx;
val(i,1) += -uvPenaltyFactor*un*ny;
if (nuH) {
double dudn = uvGrad(i,0)*nx + uvGrad(i,1)*ny;// + _UVgrad(i,2)*nz;
double dvdn = uvGrad(i,3)*nx + uvGrad(i,4)*ny;// + _UVgrad(i,5)*nz;
double nGradUVn = dudn*nx + dvdn*ny;
val(i,0) += -(*nuH)(i,0)*nGradUVn*nx;// - jumpPenaltyX;
val(i,1) += -(*nuH)(i,0)*nGradUVn*ny;// - jumpPenaltyY; // only very small effect
}
}
}
}
dgBoundaryCondition *dgConservationLawSW3dMomentum::newBoundaryWallLinear(double uvPenaltyFactor){
this->checkSetup();
return newBoundaryConditionMember(*this, &dgConservationLawSW3dMomentum::boundaryWallLinear, uvPenaltyFactor);
}
void dgConservationLawSW3dMomentum::boundarySurface(functorCache &cm, fullMatrix<double> &val, const functor *windF) const {
const fullMatrix<double> &normals = cm.get(*function::getNormals());
// const fullMatrix<double> &UV = cm.get(*function::getSolution());
......@@ -515,32 +482,6 @@ dgBoundaryCondition *dgConservationLawSW3dMomentum::newBoundarySurface(const fun
return newBoundaryConditionMember(*this, &dgConservationLawSW3dMomentum::boundarySurface, wind);
}
void dgConservationLawSW3dMomentum::boundaryOpenFree(functorCache &cm, fullMatrix<double> &val) const{
const fullMatrix<double> &normals = cm.get(*function::getNormals());
const fullMatrix<double> &eta = cm.get(*_eta,0);
const fullMatrix<double> &bath = cm.get(*_bathymetry,0);
fullMatrix<double> &externalSol = cm.set(*function::getSolution(),1);
externalSol.resize(cm.nEvaluationPoint(), 2, false);
fullMatrix<double> &externalUVInt = cm.set(*_uv2d,1);
externalUVInt.resize(cm.nEvaluationPoint(), 2, false);
cm.setDefaultToSide0(true);
double g = slim3dParameters::g;
for(int i=0;i<cm.nEvaluationPoint(); i++) {
double nx = normals(i,0);
double ny = normals(i,1);
double H = bath(i,0) + eta(i,0);
externalSol(i,0) = nx * sqrt(g/H) * eta(i,0);
externalSol(i,1) = ny * sqrt(g/H) * eta(i,0);
externalUVInt(i,0) = nx * sqrt(g/H) * eta(i,0);
externalUVInt(i,1) = ny * sqrt(g/H) * eta(i,0);
}
val.copy(cm.get(*getInterfaceTerm0("")));
}
dgBoundaryCondition *dgConservationLawSW3dMomentum::newBoundaryOpenFree(){
return newBoundaryConditionMember(*this, &dgConservationLawSW3dMomentum::boundaryOpenFree);
}
// ********* 3D horizontal momentum for U or V, vertical diffusion only
dgFaceTermAssembler *dgConservationLawSW3dMomWx::getFaceTerm(const dgGroupOfFaces &group) const
......
......@@ -10,7 +10,6 @@ class dgConservationLawSW3dMomentum : public dgConservationLawFunction {
void riemann(functorCache &cm, fullMatrix<double> &val) const;
void diffusiveFlux(functorCache &cm, fullMatrix<double> &val) const;
void maxConvectiveSpeed(functorCache &cm, fullMatrix<double> &val) const;
void boundaryWallLinear(functorCache &cm, fullMatrix<double> &val, double uvPenaltyFactor) const ;
void boundarySurface(functorCache &cm, fullMatrix<double> &val, const functor *windF) const;
const functor *_uv2d, *_uv2dGrad, *_eta, *_rGrad, *_w, *_wM, *_dwMdz, *_linearDissipation, *_quadraticDissipation, *_source, *_coriolisFactor, *_bathymetry, *_bathymetryMinCG;
const functor *_fzero, *_fzerov;
......@@ -22,7 +21,6 @@ class dgConservationLawSW3dMomentum : public dgConservationLawFunction {
void boundaryWall(functorCache &m, fullMatrix<double> &val, const boundaryWallArg arg) const;
class boundaryWallVerticalBottomArg;
void boundaryWallVerticalBottom(functorCache &m, fullMatrix<double> &val, const boundaryWallVerticalBottomArg arg) const;
void boundaryOpenFree(functorCache &cm, fullMatrix<double> &val) const;
public:
void setup();
void diffusivityTensor(functorCache &cm, fullMatrix<double> &val) const;
......@@ -75,13 +73,8 @@ class dgConservationLawSW3dMomentum : public dgConservationLawFunction {
dgBoundaryCondition *newBoundaryWall(const functor* extDiffusiveFlux=NULL, const functor* zBotDist=NULL);
/** same as boundary Wall, but for vertical bottom (uv2d is different) */
dgBoundaryCondition *newBoundaryWallVerticalBottom(const functor* extDiffusiveFlux=NULL, const functor* zBotDist=NULL);
/** Lateral wall boundary condition with normal velocit penalty
* By default uvPenaltyFactor=0 */
dgBoundaryCondition *newBoundaryWallLinear(double uvPenaltyFactor=0);
/**bnd condition for free surface */
dgBoundaryCondition *newBoundarySurface(const functor* wind);
/** open boundary condition, with etaOut = etaIn, and u = n*sqrt(gH) /H */
dgBoundaryCondition *newBoundaryOpenFree();
};
/** Vertical advection and diffusion for 3d horizontal momentum equation for u or v.
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment