XC Open source finite element analysis program
BrickUP.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 //----------------------------------------------------------------------------
27 /* ****************************************************************** **
28 ** OpenSees - Open System for Earthquake Engineering Simulation **
29 ** Pacific Earthquake Engineering Research Center **
30 ** **
31 ** **
32 ** (C) Copyright 1999, The Regents of the University of California **
33 ** All Rights Reserved. **
34 ** **
35 ** Commercial use of this program without express permission of the **
36 ** University of California, Berkeley, is strictly prohibited. See **
37 ** file 'COPYRIGHT' in main directory for information on usage and **
38 ** redistribution, and for a DISCLAIMER OF ALL WARRANTIES. **
39 ** **
40 ** Developed by: **
41 ** Frank McKenna (fmckenna@ce.berkeley.edu) **
42 ** Gregory L. Fenves (fenves@ce.berkeley.edu) **
43 ** Filip C. Filippou (filippou@ce.berkeley.edu) **
44 ** **
45 ** ****************************************************************** */
46 
48 // Description: This file contains the class declaration for //
49 // BrickUP, an 8-node cubic element for solid-fluid fully coupled analysis. //
50 // This implementation is a simplified u-p formulation of Biot theory //
51 // (u - solid displacement, p - fluid pressure). Each element node has three //
52 // DOFs for u and 1 DOF for p. //
53 // //
54 // Written by Zhaohui Yang (March 2004) //
55 // //
57 
58 // $Revision: 1.1 $
59 // $Date: 2005/09/22 21:28:36 $
60 // $Source: /usr/local/cvs/OpenSees/SRC/element/UP-ucsd/BrickUP.h,v $
61 
62 // by Zhaohui Yang (Modified based on Ed "C++" Love's Brick element)
63 //
64 // Eight node BrickUP element
65 //
66 
67 #include <domain/mesh/element/volumen/BrickBase.h>
68 #include "domain/mesh/element/utils/body_forces/BodyForces3D.h"
69 
70 namespace XC{
71 class NDMaterial;
73 //
75 class BrickUP : public BrickBase
76  {
77  private :
78  BodyForces3D bf;
79  double rho; // Fluid mass per unit volume
80  double kc; // combined bulk modulus
81  double perm[3]; // permeability
82  mutable Matrix *Ki;
83 
84  //static data
85  static Matrix stiff ;
86  static Vector resid ;
87  static Matrix mass ;
88  static Matrix damp ;
89 
90  //quadrature data
91  static const double sg[2] ;
92  static const double wg[8] ;
93 
94 
95  //local nodal coordinates, three coordinates for each of four nodes
96  // static double xl[3][8] ;
97  static double xl[][8] ;
98 
99  void formInertiaTerms(int tangFlag) const;
100  void formDampingTerms(int tangFlag) const;
101  void formResidAndTangent(int tang_flag) const;
102  double mixtureRho(int ipt) const;
103  void computeBasis(void) const;
104  const Matrix& computeBbar(int node,const double shp[4][8],const double shpBar[4][8] ) ;
105  const Matrix& computeB( int node, const double shp[4][8] ) const;
106  Matrix transpose( int dim1, int dim2, const Matrix &M ) ;
107  protected:
108  int sendData(CommParameters &);
109  int recvData(const CommParameters &);
110  public :
111  //null constructor
112  BrickUP(void);
113 
114  //full constructor
115  BrickUP( int tag, int node1,int node2,int node3,int node4,int node5,int node6,int node7,int node8,NDMaterial &theMaterial,double bulk, double rhof, double perm1, double perm2, double perm3, const BodyForces3D &bForces= BodyForces3D()) ;
116  Element *getCopy(void) const;
117  //destructor
118  virtual ~BrickUP(void) ;
119 
120  //set domain
121  void setDomain( Domain *theDomain ) ;
122 
123  //return number of dofs
124  int getNumDOF(void) const;
125 
126  //print out element data
127  void Print( std::ostream &s, int flag ) ;
128 
129  //return stiffness matrix
130  const Matrix &getTangentStiff(void) const;
131  const Matrix &getInitialStiff(void) const;
132  const Matrix &getDamp(void) const;
133  const Matrix &getMass(void) const;
134 
135  int addLoad(ElementalLoad *theLoad, double loadFactor);
136  int addInertiaLoadToUnbalance(const Vector &accel);
137 
138  //get residual
139  const Vector &getResistingForce(void) const;
140  const Vector &getResistingForceIncInertia(void) const;
141 
142  virtual int sendSelf(CommParameters &);
143  virtual int recvSelf(const CommParameters &);
144 
145  Response *setResponse(const std::vector<std::string> &argv, Information &eleInformation);
146  int getResponse(int responseID, Information &eleInformation);
147 
148 
149  // RWB; PyLiq1 & TzLiq1 need to see the excess pore pressure and initial stresses.
150  friend class PyLiq1;
151  friend class TzLiq1;
152 } ;
153 
154 } // end of XC namespace
int sendData(CommParameters &)
Send object members through the channel being passed as parameter.
Definition: BrickUP.cpp:998
??.
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
Element * getCopy(void) const
Virtual constructor.
Definition: BrickUP.cpp:127
void Print(std::ostream &s, int flag)
Imprime el objeto.
Definition: BrickUP.cpp:152
Information about an element.
Definition: Information.h:80
Hexaedro de ocho nodos.
Definition: BrickUP.h:75
Base calass for the finite elements.
Definition: Element.h:104
Base class for loads over elements.
Definition: ElementalLoad.h:73
virtual int sendSelf(CommParameters &)
Sends object through the channel being passed as parameter.
Definition: BrickUP.cpp:1018
const Vector & getResistingForceIncInertia(void) const
Returns the action of the element over its attached nodes. Computes damping matrix.
Definition: BrickUP.cpp:571
Definition: Matrix.h:82
const Matrix & getDamp(void) const
Returns the matriz de amortiguamiento.
Definition: BrickUP.cpp:380
??.
Definition: PyLiq1.h:62
Body forces over an element.
Definition: BodyForces3D.h:39
void setDomain(Domain *theDomain)
Sets the domain for the element.
Definition: BrickUP.cpp:139
virtual int recvSelf(const CommParameters &)
Receives object through the channel being passed as parameter.
Definition: BrickUP.cpp:1032
const Matrix & getMass(void) const
Returns the mass matrix.
Definition: BrickUP.cpp:369
int recvData(const CommParameters &)
Receives object members through the channel being passed as parameter.
Definition: BrickUP.cpp:1008
Communication parameters between processes.
Definition: CommParameters.h:65
Base class for hexahedra.
Definition: BrickBase.h:44
================================================================================
Definition: ContinuaReprComponent.h:34
Definition: Response.h:71