XC Open source finite element analysis program
FourNodeQuadUP.h
1 //----------------------------------------------------------------------------
2 // XC program; finite element analysis code
3 // for structural analysis and design.
4 //
5 // Copyright (C) Luis Claudio Pérez Tato
6 //
7 // This program derives from OpenSees <http://opensees.berkeley.edu>
8 // developed by the «Pacific earthquake engineering research center».
9 //
10 // Except for the restrictions that may arise from the copyright
11 // of the original program (see copyright_opensees.txt)
12 // XC is free software: you can redistribute it and/or modify
13 // it under the terms of the GNU General Public License as published by
14 // the Free Software Foundation, either version 3 of the License, or
15 // (at your option) any later version.
16 //
17 // This software is distributed in the hope that it will be useful, but
18 // WITHOUT ANY WARRANTY; without even the implied warranty of
19 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
20 // GNU General Public License for more details.
21 //
22 //
23 // You should have received a copy of the GNU General Public License
24 // along with this program.
25 // If not, see <http://www.gnu.org/licenses/>.
26 //----------------------------------------------------------------------------
28 // Description: This file contains the class declaration for FourNodeQuadUP. //
29 // FourNodeQuadUP is a 4-node plane strain element for solid-fluid fully //
30 // coupled analysis. This implementation is a simplified u-p formulation //
31 // of Biot theory (u - solid displacement, p - fluid pressure). Each element //
32 // node has two DOFs for u and 1 DOF for p. //
33 // Written by Zhaohui Yang (May 2002) //
34 // based on FourNodeQuad element by Michael Scott //
36 
37 // $Revision: 1.1 $
38 // $Date: 2005/09/22 21:28:36 $
39 // $Source: /usr/local/cvs/OpenSees/SRC/element/UP-ucsd/FourNodeQuadUP.h,v $
40 
41 #ifndef FourNodeQuadUP_h
42 #define FourNodeQuadUP_h
43 
44 #include "domain/mesh/element/plane/QuadBase4N.h"
45 #include "domain/mesh/element/utils/physical_properties/SolidMech2D.h"
46 #include "domain/mesh/element/utils/body_forces/BodyForces2D.h"
47 
48 namespace XC {
49 class Node;
50 class NDMaterial;
51 class Response;
52 
54 //
60 class FourNodeQuadUP : public QuadBase4N<SolidMech2D>
61  {
62  private:
63 
64  static Matrix K;
65  static Vector P;
66  BodyForces2D bf;
67  Vector pressureLoad;
68 
69  double rho;
70  double kc;
71  double pressure;
72  double perm[2];
73 
74  static double shp[3][4][4];
75  static double pts[4][2];
76  static double wts[4];
77  static double dvol[4];
78  static double shpBar[3][4];
79 
80  Node *nd1Ptr(void);
81  const Node *nd1Ptr(void) const;
82  Node *nd2Ptr(void);
83  const Node *nd2Ptr(void) const;
84  Node *nd3Ptr(void);
85  const Node *nd3Ptr(void) const;
86  Node *nd4Ptr(void);
87  const Node *nd4Ptr(void) const;
88 
89  // private member functions - only objects of this class can call these
90  double mixtureRho(int ipt) const; // Mixture mass density at integration point i
91  void shapeFunction(void) const;
92  void setPressureLoadAtNodes(void);
93 
94  mutable Matrix *Ki;
95  protected:
96  int sendData(CommParameters &cp);
97  int recvData(const CommParameters &cp);
98  public:
99  FourNodeQuadUP(int tag, int nd1, int nd2, int nd3, int nd4, NDMaterial &m,const std::string &type, double t, double bulk, double rhof, double perm1, double perm2,const BodyForces2D &bForces= BodyForces2D(), double p = 0.0);
100  Element *getCopy(void) const;
101  FourNodeQuadUP();
102  virtual ~FourNodeQuadUP();
103 
104  int getNumDOF(void) const;
105  void setDomain(Domain *theDomain);
106 
107  // public methods to set the state of the element
108  int update(void);
109 
110  // public methods to obtain stiffness, mass, damping and residual information
111  const Matrix &getTangentStiff(void) const;
112  const Matrix &getInitialStiff(void) const;
113  const Matrix &getDamp(void) const;
114  const Matrix &getMass(void) const;
115 
116  int addLoad(ElementalLoad *theLoad, double loadFactor);
117  int addInertiaLoadToUnbalance(const Vector &accel);
118  const Vector &getResistingForce(void) const;
119  const Vector &getResistingForceIncInertia(void) const;
120 
121  // public methods for element output
122  int sendSelf(CommParameters &);
123  int recvSelf(const CommParameters &);
124  void Print(std::ostream &s, int flag =0);
125 
126  Response *setResponse(const std::vector<std::string> &argv, Information &eleInformation);
127  int getResponse(int responseID, Information &eleInformation);
128 
129  int setParameter(const std::vector<std::string> &argv, Parameter &param);
130  int updateParameter(int parameterID, Information &info);
131 
132 
133  // RWB; PyLiq1 & TzLiq1 need to see the excess pore pressure and initial stresses.
134  friend class PyLiq1;
135  friend class TzLiq1;
136  };
137 } // end of XC namespace
138 
139 #endif
140 
void setDomain(Domain *theDomain)
Sets the domain for the element.
Definition: FourNodeQuadUP.cpp:138
??.
Definition: TzLiq1.h:60
Domain (mesh and boundary conditions) of the finite element model.
Definition: Domain.h:98
Base class for 2D and 3D materials.
Definition: NDMaterial.h:91
Definition: Vector.h:82
int recvData(const CommParameters &cp)
Receives object members through the channel being passed as parameter.
Definition: FourNodeQuadUP.cpp:614
Mesh node.
Definition: Node.h:99
const Matrix & getMass(void) const
Returns the mass matrix.
Definition: FourNodeQuadUP.cpp:355
Information about an element.
Definition: Information.h:80
Element * getCopy(void) const
Virtual constructor.
Definition: FourNodeQuadUP.cpp:110
Base calass for the finite elements.
Definition: Element.h:104
Base class for loads over elements.
Definition: ElementalLoad.h:73
const Vector & getResistingForceIncInertia(void) const
Returns the action of the element over its attached nodes. Computes damping matrix.
Definition: FourNodeQuadUP.cpp:525
Base class for 4 node quads.
Definition: QuadBase4N.h:45
int sendData(CommParameters &cp)
Send object members through the channel being passed as parameter.
Definition: FourNodeQuadUP.cpp:605
Definition: Matrix.h:82
Four-node plane-strain element using bilinear isoparametric formulation. This element is implemented ...
Definition: FourNodeQuadUP.h:60
??.
Definition: PyLiq1.h:62
Body forces over an element.
Definition: BodyForces2D.h:39
int update(void)
Actualiza el estado of the element.
Definition: FourNodeQuadUP.cpp:148
Definition: Parameter.h:65
const Matrix & getDamp(void) const
Returns the matriz de amortiguamiento.
Definition: FourNodeQuadUP.cpp:291
Communication parameters between processes.
Definition: CommParameters.h:65
================================================================================
Definition: ContinuaReprComponent.h:34
void Print(std::ostream &s, int flag=0)
Imprime el objeto.
Definition: FourNodeQuadUP.cpp:653
Definition: Response.h:71