OOFEM  2.4
OOFEM.org - Object Oriented Finite Element Solver
tr1_2d_supg2_axi.h
Go to the documentation of this file.
1 /*
2  *
3  * ##### ##### ###### ###### ### ###
4  * ## ## ## ## ## ## ## ### ##
5  * ## ## ## ## #### #### ## # ##
6  * ## ## ## ## ## ## ## ##
7  * ## ## ## ## ## ## ## ##
8  * ##### ##### ## ###### ## ##
9  *
10  *
11  * OOFEM : Object Oriented Finite Element Code
12  *
13  * Copyright (C) 1993 - 2013 Borek Patzak
14  *
15  *
16  *
17  * Czech Technical University, Faculty of Civil Engineering,
18  * Department of Structural Mechanics, 166 29 Prague, Czech Republic
19  *
20  * This library is free software; you can redistribute it and/or
21  * modify it under the terms of the GNU Lesser General Public
22  * License as published by the Free Software Foundation; either
23  * version 2.1 of the License, or (at your option) any later version.
24  *
25  * This program is distributed in the hope that it will be useful,
26  * but WITHOUT ANY WARRANTY; without even the implied warranty of
27  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
28  * Lesser General Public License for more details.
29  *
30  * You should have received a copy of the GNU Lesser General Public
31  * License along with this library; if not, write to the Free Software
32  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
33  */
34 
35 #ifndef tr1_2d_supg2_axi_h
36 #define tr1_2d_supg2_axi_h
37 
38 #include "tr1_2d_supg.h"
39 
40 #define _IFT_TR1_2D_SUPG2_AXI_Name "tr1supg2axi"
41 
42 namespace oofem {
54 {
55 protected:
60  Polygon myPoly [ 2 ];
61  std::vector< FloatArray > vcoords [ 2 ];
62 
68  int mat [ 2 ];
69 
70 public:
71  TR1_2D_SUPG2_AXI(int n, Domain * d);
72  virtual ~TR1_2D_SUPG2_AXI();
73 
74  virtual void computeAccelerationTerm_MB(FloatMatrix &answer, TimeStep *tStep);
75  virtual void computeAdvectionTerm_MB(FloatArray &answer, TimeStep *tStep);
76  virtual void computeAdvectionDerivativeTerm_MB(FloatMatrix &answer, TimeStep *tStep);
77  virtual void computeDiffusionTerm_MB(FloatArray &answer, TimeStep *tStep);
78  virtual void computeDiffusionDerivativeTerm_MB(FloatMatrix &answer, MatResponseMode mode, TimeStep *tStep);
79  virtual void computePressureTerm_MB(FloatMatrix &answer, TimeStep *tStep);
80  virtual void computeLSICStabilizationTerm_MB(FloatMatrix &answer, TimeStep *tStep);
81  virtual void computeLinearAdvectionTerm_MC(FloatMatrix &answer, TimeStep *tStep);
82  virtual void computeAdvectionTerm_MC(FloatArray &answer, TimeStep *tStep);
83  virtual void computeAdvectionDerivativeTerm_MC(FloatMatrix &answer, TimeStep *tStep);
84  virtual void computeDiffusionDerivativeTerm_MC(FloatMatrix &answer, TimeStep *tStep);
85  virtual void computeDiffusionTerm_MC(FloatArray &answer, TimeStep *tStep);
86  virtual void computeAccelerationTerm_MC(FloatMatrix &answer, TimeStep *tStep);
87  virtual void computePressureTerm_MC(FloatMatrix &answer, TimeStep *tStep);
88  virtual void computeBCRhsTerm_MB(FloatArray &answer, TimeStep *tStep);
89  virtual void computeBCRhsTerm_MC(FloatArray &answer, TimeStep *tStep);
90 
91  virtual void updateStabilizationCoeffs(TimeStep *tStep);
93  virtual double computeCriticalTimeStep(TimeStep *tStep);
94 
95  // definition
96  virtual const char *giveClassName() const { return "TR1_2D_SUPG2_AXI"; }
97  virtual const char *giveInputRecordName() const { return _IFT_TR1_2D_SUPG2_AXI_Name; }
98  virtual MaterialMode giveMaterialMode() { return _2dAxiFlow; }
100  virtual void giveInputRecord(DynamicInputRecord &input);
101 
102  virtual void printOutputAt(FILE *file, TimeStep *tStep);
103 
104 #ifdef __OOFEG
106  int node, TimeStep *tStep);
107  // Graphics output
108  virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep);
109  virtual void drawScalar(oofegGraphicContext &gc, TimeStep *tStep);
110  //virtual void drawDeformedGeometry(oofegGraphicContext &gc, TimeStep *tStep, UnknownType) {}
111 #endif
112 
113 protected:
114  virtual void computeGaussPoints();
115  virtual void computeDeviatoricStress(FloatArray &answer, GaussPoint *gp, TimeStep *);
116  void updateVolumePolygons(Polygon &referenceFluidPoly, Polygon &secondFluidPoly, int &rfPoints, int &sfPoints,
117  const FloatArray &normal, const double p, bool updFlag);
118  double computeVolumeAroundID(GaussPoint *gp, integrationDomain id, const std::vector< FloatArray > &idpoly);
119  double computeRadiusAt(GaussPoint *gp);
120  void computeBMtrx(FloatMatrix &answer, GaussPoint *gp);
121  void computeNVector(FloatArray &answer, GaussPoint *gp);
122  void updateIntegrationRules();
123  Material *_giveMaterial(int indx) { return domain->giveMaterial(mat [ indx ]); }
124 
125  virtual void NodalAveragingRecoveryMI_computeNodalValue(FloatArray &answer, int node,
126  InternalStateType type, TimeStep *tStep);
127 
129  virtual void SPRNodalRecoveryMI_giveDofMansDeterminedByPatch(IntArray &answer, int pap);
130  virtual int SPRNodalRecoveryMI_giveNumberOfIP();
132 
133  virtual double computeLEPLICVolumeFraction(const FloatArray &n, const double p, LEPlic *matInterface, bool updFlag);
134  virtual void formMaterialVolumePoly(Polygon &matvolpoly, LEPlic *matInterface,
135  const FloatArray &normal, const double p, bool updFlag);
136  virtual void formVolumeInterfacePoly(Polygon &matvolpoly, LEPlic *matInterface,
137  const FloatArray &normal, const double p, bool updFlag);
138  virtual double truncateMatVolume(const Polygon &matvolpoly, double &volume);
139  virtual void giveElementCenter(LEPlic *mat_interface, FloatArray &center, bool updFlag);
140  virtual void formMyVolumePoly(Polygon &myPoly, LEPlic *mat_interface, bool updFlag);
141  virtual Element *giveElement() { return this; }
142  virtual double computeMyVolume(LEPlic *matInterface, bool updFlag);
143 };
144 } // end namespace oofem
145 #endif // tr1_2d_supg2_axi_h
InternalStateType
Type representing the physical meaning of element or constitutive model internal variable.
integrationDomain
Used by integrator class to supply integration points for proper domain to be integrated (Area...
virtual IRResultType initializeFrom(InputRecord *ir)
Initializes receiver according to object description stored in input record.
Class and object Domain.
Definition: domain.h:115
virtual void computeDiffusionDerivativeTerm_MC(FloatMatrix &answer, TimeStep *tStep)
Computes diffusion derivative terms for mass conservation equation.
Domain * domain
Link to domain object, useful for communicating with other FEM components.
Definition: femcmpnn.h:82
int giveInternalStateAtNode(FloatArray &answer, InternalStateType type, InternalStateMode mode, int node, TimeStep *tStep)
Returns internal state variable (like stress,strain) at node of element in Reduced form...
virtual void computeLSICStabilizationTerm_MB(FloatMatrix &answer, TimeStep *tStep)
Computes SLIC stabilization term for momentum balance equation(s).
virtual void computeAdvectionDerivativeTerm_MB(FloatMatrix &answer, TimeStep *tStep)
Computes the derivative of advection terms for momentum balance equations(s) with respect to nodal ve...
virtual void formMyVolumePoly(Polygon &myPoly, LEPlic *mat_interface, bool updFlag)
Assembles receiver volume.
virtual void printOutputAt(FILE *file, TimeStep *tStep)
Prints output of receiver to stream, for given time step.
virtual void updateElementForNewInterfacePosition(TimeStep *tStep)
oofem::oofegGraphicContext gc[OOFEG_LAST_LAYER]
virtual void computeDiffusionTerm_MC(FloatArray &answer, TimeStep *tStep)
Computes diffusion terms for mass conservation equation.
virtual void NodalAveragingRecoveryMI_computeNodalValue(FloatArray &answer, int node, InternalStateType type, TimeStep *tStep)
Computes the element value in given node.
virtual void SPRNodalRecoveryMI_giveDofMansDeterminedByPatch(IntArray &answer, int pap)
Abstract base class for all finite elements.
Definition: element.h:145
void computeNVector(FloatArray &answer, GaussPoint *gp)
virtual void computeGaussPoints()
Initializes the array of integration rules member variable.
std::vector< FloatArray > vcoords[2]
virtual void computePressureTerm_MB(FloatMatrix &answer, TimeStep *tStep)
Computes pressure terms for momentum balance equations(s).
MaterialMode
Type representing material mode of integration point.
Definition: materialmode.h:89
Class implementing an array of integers.
Definition: intarray.h:61
MatResponseMode
Describes the character of characteristic material matrix.
void updateVolumePolygons(Polygon &referenceFluidPoly, Polygon &secondFluidPoly, int &rfPoints, int &sfPoints, const FloatArray &normal, const double p, bool updFlag)
virtual void drawRawGeometry(oofegGraphicContext &gc, TimeStep *tStep)
virtual void formMaterialVolumePoly(Polygon &matvolpoly, LEPlic *matInterface, const FloatArray &normal, const double p, bool updFlag)
Assembles the true element material polygon (takes receiver vof into accout).
virtual const char * giveClassName() const
virtual void computeAccelerationTerm_MC(FloatMatrix &answer, TimeStep *tStep)
Computes acceleration terms for mass conservation equation.
virtual void updateStabilizationCoeffs(TimeStep *tStep)
Updates the stabilization coefficients used for CBS and SUPG algorithms.
Material * giveMaterial(int n)
Service for accessing particular domain material model.
Definition: domain.C:281
virtual void computeAdvectionTerm_MC(FloatArray &answer, TimeStep *tStep)
Computes advection terms for mass conservation equation.
virtual SPRPatchType SPRNodalRecoveryMI_givePatchType()
virtual void computeAdvectionDerivativeTerm_MC(FloatMatrix &answer, TimeStep *tStep)
Computes the derivative of advection terms for mass conservation equation with respect to nodal veloc...
virtual MaterialMode giveMaterialMode()
Returns material mode for receiver integration points.
virtual void computeLinearAdvectionTerm_MC(FloatMatrix &answer, TimeStep *tStep)
Computes the linear advection term for mass conservation equation.
virtual void computeDiffusionTerm_MB(FloatArray &answer, TimeStep *tStep)
Computes diffusion terms for momentum balance equations(s).
Abstract base class for all material models.
Definition: material.h:95
virtual void computeBCRhsTerm_MC(FloatArray &answer, TimeStep *tStep)
Computes Rhs terms due to boundary conditions.
#define _IFT_TR1_2D_SUPG2_AXI_Name
virtual double truncateMatVolume(const Polygon &matvolpoly, double &volume)
Truncates given material polygon to receiver.
virtual double computeCriticalTimeStep(TimeStep *tStep)
Computes the critical time increment.
virtual void computeBCRhsTerm_MB(FloatArray &answer, TimeStep *tStep)
Computes Rhs terms due to boundary conditions.
virtual void SPRNodalRecoveryMI_giveSPRAssemblyPoints(IntArray &pap)
virtual double computeMyVolume(LEPlic *matInterface, bool updFlag)
Computes the volume of receiver.
virtual void computeDeviatoricStress(FloatArray &answer, GaussPoint *gp, TimeStep *)
double computeRadiusAt(GaussPoint *gp)
Abstract base class representing Lagrangian-Eulerian (moving) material interfaces.
Definition: leplic.h:153
virtual void computeDiffusionDerivativeTerm_MB(FloatMatrix &answer, MatResponseMode mode, TimeStep *tStep)
Computes the derivative of diffusion terms for momentum balance equations(s) with respect to nodal ve...
Class representing vector of real numbers.
Definition: floatarray.h:82
void computeBMtrx(FloatMatrix &answer, GaussPoint *gp)
virtual double computeLEPLICVolumeFraction(const FloatArray &n, const double p, LEPlic *matInterface, bool updFlag)
Computes corresponding volume fraction to given interface position.
Class representing 2D polygon.
Definition: geotoolbox.h:91
Implementation of matrix containing floating point numbers.
Definition: floatmatrix.h:94
IRResultType
Type defining the return values of InputRecord reading operations.
Definition: irresulttype.h:47
Class representing the general Input Record.
Definition: inputrecord.h:101
Class representing 2d linear axisymmetric triangular element for solving incompressible fluid with SU...
virtual void computeAdvectionTerm_MB(FloatArray &answer, TimeStep *tStep)
Computes nonlinear advection terms for momentum balance equations(s).
virtual void formVolumeInterfacePoly(Polygon &matvolpoly, LEPlic *matInterface, const FloatArray &normal, const double p, bool updFlag)
Assembles receiver material polygon based solely on given interface line.
Polygon myPoly[2]
myPoly[0] Occupied by reference fluid.
virtual int SPRNodalRecoveryMI_giveNumberOfIP()
virtual void computeAccelerationTerm_MB(FloatMatrix &answer, TimeStep *tStep)
Computes acceleration terms (generalized mass matrix with stabilization terms) for momentum balance e...
Class representing the a dynamic Input Record.
int mat[2]
mat[0] reference fluid mat[1] second fluid
Material * _giveMaterial(int indx)
virtual void giveInputRecord(DynamicInputRecord &input)
Setups the input record string of receiver.
Class representing 2d linear triangular element for solving incompressible fluid with SUPG solver...
Definition: tr1_2d_supg.h:68
virtual const char * giveInputRecordName() const
double p
Line constant of line segment representing interface.
Definition: leplic.h:65
TR1_2D_SUPG2_AXI(int n, Domain *d)
double computeVolumeAroundID(GaussPoint *gp, integrationDomain id, const std::vector< FloatArray > &idpoly)
the oofem namespace is to define a context or scope in which all oofem names are defined.
virtual void giveElementCenter(LEPlic *mat_interface, FloatArray &center, bool updFlag)
Computes the receiver center (in updated Lagrangian configuration).
virtual void drawScalar(oofegGraphicContext &gc, TimeStep *tStep)
Class representing integration point in finite element program.
Definition: gausspoint.h:93
FloatArray normal
Interface segment normal.
Definition: leplic.h:67
Class representing solution step.
Definition: timestep.h:80
virtual Element * giveElement()
Return number of receiver&#39;s element.
InternalStateMode
Determines the mode of internal variable.
virtual void computePressureTerm_MC(FloatMatrix &answer, TimeStep *tStep)
Computes pressure terms for mass conservation equation.

This page is part of the OOFEM documentation. Copyright (c) 2011 Borek Patzak
Project e-mail: info@oofem.org
Generated at Tue Jan 2 2018 20:07:32 for OOFEM by doxygen 1.8.11 written by Dimitri van Heesch, © 1997-2011