#include <quadtet.h>
Public Member Functions | |
double | approx (double xi, double eta, double zeta, vector &nodval) |
computes approximated quantity given by the nodal values at the given point | |
void | aver_strains (long lcid, long eid, long ri, long ci, vector &averstra, double &volume) |
computes averaged strains | |
void | bf_matrix (matrix &n, double xi, double eta, double zeta) |
computes base functions matrix | |
void | compute_eigstress (long lcid, long eid, long ri, long ci) |
computes stresses due to eigenstrains | |
void | compute_nlstress (long lcid, long eid, long ri, long ci) |
computes stresses in the integration points, laws of specified material models is used for stress evaluation | |
void | compute_nlstressincr (long lcid, long eid, long ri, long ci) |
computes increments of stresses in the integration points, specified material model is used for stress evaluation | |
void | compute_nonloc_nlstress (long lcid, long eid, long ri, long ci) |
computes stresses in the integration points, nonlocal material model is used for stress evaluation | |
void | eigstrain_forces (long lcid, long eid, long ri, long ci, vector &nfor, vector &x, vector &y, vector &z) |
computes nodal forces of the element caused by eigenstrains | |
void | elem_integration (integratedquant iq, long lcid, long eid, long ri, long ci, vector &nv, vector &x, vector &y, vector &z) |
integrates required quantity over the element | |
void | geom_matrix (matrix &gm, vector &x, vector &y, vector &z, double xi, double eta, double zeta, double &jac) |
computes geometrical matrix | |
void | geom_matrix_block (matrix &gm, vector &x, vector &y, vector &z, double xi, double eta, double zeta, double &jac) |
computes block of geometrical matrix | |
void | incr_internal_forces (long lcid, long eid, long ri, long ci, vector &ifor, vector &x, vector &y, vector &z) |
computes increments of internal forces of the element | |
void | inicipval (long eid, long ri, long ci, matrix &nodval, inictype *ictn) |
sets internal data of integration points to specified initial values | |
void | internal_forces (long lcid, long eid, long ri, long ci, vector &ifor, vector &x, vector &y, vector &z) |
computes internal forces of the element | |
void | intpointval (long eid, vector &nodval, vector &ipval) |
approximates nodal values to the the integration point values (quadratic approximation is used) | |
void | intpointval2 (long eid, vector &nodval, vector &ipval) |
approximates nodal values to the the integration point values (linear approximation is used) | |
void | ip_elast_stresses (long lcid, long eid, long ri, long ci) |
computes stresses in the integration points, elastic law is used | |
void | ip_strains (long lcid, long eid, long ri, long ci, vector &x, vector &y, vector &z, vector &r) |
computes strains in the integration points | |
void | ip_stresses (long lcid, long eid, long ri, long ci) |
computes stresses in the integration points | |
void | ipcoord (long eid, long ipp, long ri, long ci, vector &coord) |
returns natural coordinates of integration points of the element | |
void | ipvolume (long eid, long ri, long ci) |
computes volumes appropriated to integration points | |
void | load_matrix (long eid, matrix &lm) |
computes load matrix | |
void | local_values (long lcid, long eid, long ri, long ci) |
computes local values of averaged quantities (used in nonlocal material models) | |
void | locglob_nodeval (long is, vector &nv, double *tnv, vector &x, vector &y, vector &z) |
transforms surface load given in the local coordinate system to the global coordinate system | |
void | mass_matrix (long eid, matrix &mm) |
computes mass matrix | |
void | nod_eqother_ip (long eid, long ri, long ci) |
copies other values from integration point to nodes | |
void | nod_strains_comp (long lcid, long eid, double **stra) |
computes strains at nodes | |
void | nod_strains_ip (long lcid, long eid, long ri, long ci) |
copies strains from integration point to nodes | |
void | nod_stresses_comp (long lcid, long eid, long ri, long ci, double **stra, double **stre) |
computes stresses at nodes | |
void | nod_stresses_ip (long lcid, long eid, long ri, long ci) |
copies stresses from integration point to nodes | |
void | node_forces_surf (long lcid, long eid, long *is, double *nv, vector &nf) |
computes nodal forces caused by surface load | |
void | nonloc_internal_forces (long lcid, long eid, long ri, long ci, vector &ifor, vector &x, vector &y, vector &z) |
computes internal forces of the element nonlocally | |
quadtet (void) | |
void | res_eigstrain_forces (long lcid, long eid, vector &nfor) |
computes resulting nodal forces of the element caused by eigenstrains | |
void | res_incr_internal_forces (long lcid, long eid, vector &ifor) |
computes resulting increments of internal forces of the element | |
void | res_internal_forces (long lcid, long eid, vector &ifor) |
computes resulting internal forces of the element | |
void | res_ip_strains (long lcid, long eid) |
computes resulting strains in the integration points | |
void | res_ip_stresses (long lcid, long eid) |
computes resulting stresses in the integration points | |
void | res_load_matrix (long eid, matrix &lm) |
computes resulting load matrix | |
void | res_mass_matrix (long eid, matrix &mm) |
computes resulting mass matrix | |
void | res_nonloc_internal_forces (long lcid, long eid, vector &ifor) |
computes resulting internal forces of the element nonlocally | |
void | res_stiffness_matrix (long eid, matrix &sm) |
computes resulting stiffness matrix | |
void | stiffness_matrix (long eid, long ri, long ci, matrix &sm) |
computes stiffness matrix | |
void | strains (long lcid, long eid, long ri, long ci) |
computes required strains | |
void | stresses (long lcid, long eid, long ri, long ci) |
computes required stresses | |
void | transf_matrix (ivector &nodes, matrix &tmat) |
computes transformation matrix | |
double | volumeip (long eid, double w) |
computes volume appropriate to integration point | |
~quadtet (void) | |
Public Attributes | |
long * | cncomp |
array containing cumulative numbers of components of stress and strain tensors | |
long | intordb |
array of numbers of integration points in sets | |
long | intordmm |
order of integration of mass matrix | |
long ** | intordsm |
array of orders of integration of stiffness matrix | |
long | napfun |
number of approximated functions on the element | |
long | nb |
number of blocks | |
long * | ncomp |
array containing numbers of components of stress and strain tensors | |
long | ndofe |
number of DOFs on the element | |
long | ned |
number of edges | |
long ** | nip |
array of numbers of integration points on surface | |
long | nne |
number of nodes on one element | |
long | nned |
number of nodes on one edge | |
long | nnsurf |
number of nodes on one surface | |
long | nsurf |
number of surfaces | |
strastrestate | ssst |
stress/strain state | |
long | tncomp |
number of components of the strain and stress tensors | |
long | tnip |
total number of integration points on element |
Class quadtet defines tetraheadronal finite element with tri-quadratic approximation functions.
Created by Tomas Koudelka
Definition at line 16 of file quadtet.h.
quadtet | ( | void | ) |
~quadtet | ( | void | ) |
double approx | ( | double | xi, | |
double | eta, | |||
double | zeta, | |||
vector & | nodval | |||
) |
computes approximated quantity given by the nodal values at the given point
Function approximates quantity defined by nodal values to the given point.
xi,eta,zeta | - natural coordinates of required point | |
nodval | - nodal values |
Created by Tomas Koudelka, 11.2008
Definition at line 119 of file quadtet.cpp.
References vector::a, bf_quad_tet(), f, nne, and scprd().
Referenced by inicipval(), interpolate(), intpointval(), ipcoord(), and mass_matrix().
void aver_strains | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci, | |||
vector & | averstra, | |||
double & | volume | |||
) |
computes averaged strains
Functiion computes averaged strains at the given element.
lcid | - load case id | |
eid | - element id | |
ri | - row index of the given integration point block | |
ci | - column index of the given integration point block | |
averstra | - vector of avreaged values of strains (have to be divided by element volume) (output) | |
volume | - element volume (output) |
Created by Tomas Koudelka, 1.2009
Definition at line 2150 of file quadtet.cpp.
References vector::a, allocv(), destrv(), mechtop::elements, gauss_points_tet(), mechtop::give_node_coord3d(), mechmat::givestrain(), intordsm, element::ipp, jac_3d(), Mm, Mt, vector::n, nb, ncomp, and nne.
void bf_matrix | ( | matrix & | n, | |
double | xi, | |||
double | eta, | |||
double | zeta | |||
) |
computes base functions matrix
Function assembles matrix of base functions in the required point.
n | - matrix of base functions (output) | |
xi,eta,zeta | - natural coordinates of required point |
Created by Tomas Koudelka, 11.2008
Definition at line 143 of file quadtet.cpp.
References vector::a, bf_quad_tet(), fillm(), and nne.
Referenced by load_matrix(), mass_matrix(), and node_forces_surf().
void compute_eigstress | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes stresses due to eigenstrains
Function computes nonlocal correct stresses at integration points on element.
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 1487 of file quadtet.cpp.
References mechtop::elements, mechmat::giveeigstrain(), intordsm, element::ipp, mechmat::matstiff(), Mm, Mt, mxv(), mechmat::storeeigstress(), and tncomp.
Referenced by eigstrain_forces().
void compute_nlstress | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes stresses in the integration points, laws of specified material models is used for stress evaluation
Function computes correct stresses at integration points on element.
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 1374 of file quadtet.cpp.
References mechmat::computenlstresses(), mechtop::elements, intordsm, element::ipp, Mm, Mp, Mt, and probdesc::strcomp.
Referenced by elem_nlstresses(), and internal_forces().
void compute_nlstressincr | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes increments of stresses in the integration points, specified material model is used for stress evaluation
Function computes correct increments of stresses at integration points on element.
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 1402 of file quadtet.cpp.
References mechmat::computenlstressesincr(), mechtop::elements, intordsm, element::ipp, Mm, Mp, Mt, and probdesc::strcomp.
Referenced by incr_internal_forces().
void compute_nonloc_nlstress | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes stresses in the integration points, nonlocal material model is used for stress evaluation
Function computes nonlocal correct stresses at integration points on element
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 1459 of file quadtet.cpp.
References mechmat::compnonloc_nlstresses(), mechtop::elements, intordsm, element::ipp, Mm, Mp, Mt, and probdesc::strcomp.
Referenced by nonloc_internal_forces().
void eigstrain_forces | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci, | |||
vector & | nfor, | |||
vector & | x, | |||
vector & | y, | |||
vector & | z | |||
) |
computes nodal forces of the element caused by eigenstrains
Function computes contributions from eigenstrains to the right hand side
lcid | - load case id | |
eid | - element id | |
ri,ci | - row and column indices | |
ifor | - vector of nodal forces due to eigenstrains | |
x,y,z | - vectors of nodal coordinates |
Created by Tomas Koudelka, 11.2008
Definition at line 1204 of file quadtet.cpp.
References compute_eigstress(), probdesc::eigstrcomp, eigstress, elem_integration(), and Mp.
Referenced by res_eigstrain_forces().
void elem_integration | ( | integratedquant | iq, | |
long | lcid, | |||
long | eid, | |||
long | ri, | |||
long | ci, | |||
vector & | nv, | |||
vector & | x, | |||
vector & | y, | |||
vector & | z | |||
) |
integrates required quantity over the element
Function integrates selected quantity over the finite element.
iq | - type of integrated quantity (see alias.h) | |
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices | |
nv | - nodal values | |
x,y,z | - node coordinates |
Created by Tomas Koudelka, 11.2008
Definition at line 1522 of file quadtet.cpp.
References vector::a, addv(), allocv(), cmulv(), cncomp, destrv(), mechtop::elements, fillv(), gauss_points_tet(), geom_matrix(), mechmat::givequantity(), intordsm, element::ipp, Mm, Mt, mtxv(), ndofe, and tncomp.
Referenced by eigstrain_forces(), incr_internal_forces(), internal_forces(), and nonloc_internal_forces().
void geom_matrix | ( | matrix & | gm, | |
vector & | x, | |||
vector & | y, | |||
vector & | z, | |||
double | xi, | |||
double | eta, | |||
double | zeta, | |||
double & | jac | |||
) |
computes geometrical matrix
Function computes strain-displacement (geometrical) matrix in the given point.
gm | - geometric matrix | |
x,y,z | - vectors containing element node coordinates | |
xi,eta,zeta | - natural coordinates of required point | |
jac | - Jacobian (output) |
Created by Tomas Koudelka, 11.2008
Definition at line 174 of file quadtet.cpp.
References vector::a, derivatives_3d(), dx_bf_quad_tet(), dy_bf_quad_tet(), dz_bf_quad_tet(), fillm(), and nne.
Referenced by elem_integration(), ip_strains(), nod_strains_comp(), and stiffness_matrix().
void geom_matrix_block | ( | matrix & | gm, | |
vector & | x, | |||
vector & | y, | |||
vector & | z, | |||
double | xi, | |||
double | eta, | |||
double | zeta, | |||
double & | jac | |||
) |
computes block of geometrical matrix
void incr_internal_forces | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci, | |||
vector & | ifor, | |||
vector & | x, | |||
vector & | y, | |||
vector & | z | |||
) |
computes increments of internal forces of the element
Function computes increments of internal forces.
lcid | - load case id | |
eid | - element id | |
ri,ci | - row and column indices | |
ifor | - vector of increments of internal forces (output) | |
x,y,z | - vectors of nodal coordinates |
Created by Tomas Koudelka, 11.2008
Definition at line 1178 of file quadtet.cpp.
References compute_nlstressincr(), elem_integration(), and stressincr.
Referenced by res_incr_internal_forces().
sets internal data of integration points to specified initial values
The function computes initial values of the given quantities at each integration point of the element from the nodal values given by the parameter nodval. Initial condition types must be the same for all nodes of the element.
eid | - element id | |
ri | - block row index | |
ci | - block column index | |
nodval | - nodal values of particular initial conditions. nodval[i][j] represents value of j-th initial condition at i-th node of the given element. | |
ictn | - array of types of initial condition for each node of element. The type of initial condition determines which values are being specified in the node. (ictn[i] & inistrain) returns nonzero if nodal values of initial strains are specified (ictn[i] & inistress) returns nonzero if nodal values of initial stresses are specified (ictn[i] & iniother) returns nonzero if nodal values of initial values of eqother array are specified (ictn[i] & inicond) returns nonzero if nodal values of other initial conditions are specified |
The function does not return anything.
Created by Tomas Koudelka 2004 Revised by Tomas Koudelka 03.2012
Definition at line 1624 of file quadtet.cpp.
References vector::a, allocv(), approx(), destrv(), mechtop::elements, intpoints::eqother, gauss_points_tet(), mechmat::ic, inicond, iniother, inistrain, inistress, intordsm, mechmat::ip, element::ipp, memset(), Mm, Mt, matrix::n, nb, intpoints::ncompeqother, intpoints::ncompstr, nne, print_err(), intpoints::strain, and intpoints::stress.
Referenced by mechbclc::inicipval().
void internal_forces | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci, | |||
vector & | ifor, | |||
vector & | x, | |||
vector & | y, | |||
vector & | z | |||
) |
computes internal forces of the element
Function computes internal forces in the case of geometrical linear computation.
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices | |
ifor | - vector of internal forces (output) | |
x,y,z | - vectors of nodal coordinates |
Created by Tomas Koudelka, 11.2008
Definition at line 1123 of file quadtet.cpp.
References compute_nlstress(), elem_integration(), and locstress.
Referenced by res_internal_forces().
approximates nodal values to the the integration point values (quadratic approximation is used)
Function interpolates the nodal values to the integration points on the element quadratic approximation functions are used
eid | - element id | |
nodval | - nodal values | |
ipval | - value at integration points |
Created by Tomas Koudelka, 11.2008
Definition at line 2042 of file quadtet.cpp.
References vector::a, allocv(), approx(), destrv(), gauss_points_tet(), intordsm, and nb.
Referenced by elem_intpointval().
approximates nodal values to the the integration point values (linear approximation is used)
Function interpolates the nodal values to the integration points on the element linear approximation functions are used
eid | - element id | |
nodval | - nodal values | |
ipval | - value at integration points |
Created by Tomas Koudelka, 11.2008
Definition at line 2092 of file quadtet.cpp.
References vector::a, allocv(), lintet::approx_nat(), destrv(), gauss_points_tet(), intordsm, Ltet, nb, and lintet::nne.
Referenced by elem_intpointval2().
void ip_elast_stresses | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes stresses in the integration points, elastic law is used
Function computes stresses at integration points of element stresses are computed from strains with the help of elastic stiffness.
lcid | - load case id | |
eid | - element id | |
ri | - row index | |
ci | - column index |
Created by Tomas Koudelka, 11.2008
Definition at line 860 of file quadtet.cpp.
References mechtop::elements, mechmat::givestrain(), intordsm, element::ipp, mechmat::matstiff(), Mm, Mt, mxv(), ndofe, mechmat::storestress(), and tncomp.
void ip_strains | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci, | |||
vector & | x, | |||
vector & | y, | |||
vector & | z, | |||
vector & | r | |||
) |
computes strains in the integration points
Function computes strains at integration points of given element. Used displacements are passed over parameter r.
lcid | - load case id | |
eid | - element id | |
ri | - row index | |
ci | - column index | |
x,y,z | - vectors of nodal coordinates | |
r | - vector of nodal displacements |
Created by Tomas Koudelka, 11.2008
Definition at line 555 of file quadtet.cpp.
References vector::a, allocm(), allocv(), cncomp, destrm(), destrv(), mechtop::elements, gauss_points_tet(), geom_matrix(), intordsm, element::ipp, Mm, Mt, mxv(), nb, ncomp, ndofe, nne, nodes, and mechmat::storestrain().
Referenced by res_ip_strains().
void ip_stresses | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes stresses in the integration points
Function computes stresses at integration points of element stresses are computed by material models.
lcid | - load case id | |
eid | - element id | |
ri | - row index | |
ci | - column index |
Created by Tomas Koudelka, 11.2008
Definition at line 830 of file quadtet.cpp.
References mechmat::computenlstresses(), mechtop::elements, intordsm, element::ipp, Mm, Mp, Mt, and probdesc::strcomp.
Referenced by res_ip_stresses().
void ipcoord | ( | long | eid, | |
long | ipp, | |||
long | ri, | |||
long | ci, | |||
vector & | coord | |||
) |
returns natural coordinates of integration points of the element
Function returns coordinates of integration points
eid | - element id | |
ipp | - integration point pointer | |
ri | - row index | |
ci | - column index | |
coord | - vector of coordinates (output) |
Function | returns coordinates of integration points in the parameter coord. |
Created by Tomas Koudelka, 11.2008
Definition at line 1574 of file quadtet.cpp.
References vector::a, approx(), mechtop::elements, gauss_points_tet(), mechtop::give_node_coord3d(), intordsm, element::ipp, Mt, and nne.
Referenced by ipcoord().
void ipvolume | ( | long | eid, | |
long | ri, | |||
long | ci | |||
) |
computes volumes appropriated to integration points
Function computes volume appropriate to integration point.
eid | - element id | |
ri | - row index of the given integration point block | |
ci | - column index of the given integration point block |
Created by Tomas Koudelka, 11.2008
Definition at line 1749 of file quadtet.cpp.
References vector::a, allocv(), destrv(), mechtop::elements, gauss_points_tet(), mechtop::give_node_coord3d(), intordsm, element::ipp, jac_3d(), Mm, Mt, nb, nne, and mechmat::storeipvol().
Referenced by ipvolume().
void load_matrix | ( | long | eid, | |
matrix & | lm | |||
) |
computes load matrix
Function computes load matrix without possible transformations at nodes.
eid | - number of element | |
lm | - load matrix |
Created by Tomas Koudelka, 11.2008
Definition at line 420 of file quadtet.cpp.
References matrix::a, vector::a, bf_matrix(), fillm(), gauss_points_tet(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), intordmm, jac_3d(), matrix::m, Mt, matrix::n, napfun, ndofe, nne, nnj(), and nodes.
Referenced by loadmat(), and res_load_matrix().
void local_values | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes local values of averaged quantities (used in nonlocal material models)
Function computes local values which will be used for averageing in nonlocal models at integration points. Mp->nonlocphase have to be 1.
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 1431 of file quadtet.cpp.
References mechmat::computenlstresses(), mechtop::elements, intordsm, element::ipp, Mm, Mp, Mt, and probdesc::strcomp.
Referenced by elem_local_values().
transforms surface load given in the local coordinate system to the global coordinate system
The function transforms nodal values of surface load expressed in local coordinate system of the given surface is to the global coordinate system.
is | - index of surface | |
nv | - vector of nodal values of surface load (3 components of load are assumed at ech node) | |
tnv | - resulting surface load in the globla coordinate system (output) | |
x | - nodal x-coordinates | |
y | - nodal y-coordinates | |
z | - nodal z-coordinates |
Created by Tomas Koudelka,
Definition at line 1947 of file quadtet.cpp.
References vector::a, matrix::a, cmulv(), mxv(), and scprd().
Referenced by node_forces_surf().
void mass_matrix | ( | long | eid, | |
matrix & | mm | |||
) |
computes mass matrix
Function computes mass matrix. If it is required, nodal values are transformed to the local coordinate systems.
eid | - number of element | |
mm | - mass matrix |
Created by Tomas Koudelka, 11.2008
Definition at line 347 of file quadtet.cpp.
References matrix::a, vector::a, approx(), bf_matrix(), fillm(), gauss_points_tet(), mechcrsec::give_density(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), intordmm, jac_3d(), matrix::m, Mc, Mt, matrix::n, napfun, ndofe, nne, nnj(), and nodes.
Referenced by massmat(), and res_mass_matrix().
void nod_eqother_ip | ( | long | eid, | |
long | ri, | |||
long | ci | |||
) |
copies other values from integration point to nodes
Function computes other values in nodes of element.
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 1052 of file quadtet.cpp.
References vector::a, allocv(), cmulv(), destrv(), mechtop::elements, gauss_points_tet(), mechtop::give_elemnodes(), mechmat::giveeqother(), mechmat::givencompeqother(), intordsm, element::ipp, Mm, Mp, Mt, nne, nod, mechtop::nodes, nodip_quadtet(), probdesc::otheraver, node::storeother(), and volumeip().
Referenced by compute_nodeothers().
void nod_strains_comp | ( | long | lcid, | |
long | eid, | |||
double ** | stra | |||
) |
computes strains at nodes
Function computes strains at nodes of element.
lcid | - load case id | |
eid | - element id | |
stra | - array for strain components (output) |
stra[i][j] - the j-th strain component at the i-th node
Created by Tomas Koudelka, 11.2008
Definition at line 671 of file quadtet.cpp.
References allocm(), allocv(), copyv(), destrm(), destrv(), eldispl(), geom_matrix(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), lgvectortransf(), mechtop::locsystems(), Mt, mxv(), ndofe, nne, nodcoord_quadtet(), nodes, tncomp, and transf_matrix().
Referenced by strains().
void nod_strains_ip | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
copies strains from integration point to nodes
Function computes strains at nodes of given element. The strains are copied to the node from the nearest integration point and they are averaged if it is prescribed.
lcid | - load case id | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 609 of file quadtet.cpp.
References vector::a, allocv(), cmulv(), destrv(), mechtop::elements, gauss_points_tet(), mechtop::give_elemnodes(), mechmat::givestrain(), intordsm, element::ipp, Mm, Mp, Mt, nne, nod, mechtop::nodes, nodip_quadtet(), node::storestrain(), probdesc::strainaver, tncomp, and volumeip().
Referenced by compute_nodestrains(), and strains().
void nod_stresses_comp | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci, | |||
double ** | stra, | |||
double ** | stre | |||
) |
computes stresses at nodes
Function computes stresses in nodes.
lcid | - load case id | |
eid | - element id | |
ri,ci | - row and column indices | |
stra | - strains at nodes, stra[i][j] means eps[j] at node i | |
stre | - stresses at nodes - output parameter, stre[i][j] means sig[j] at node i (output) |
Created by Tomas Koudelka, 11.2008
Definition at line 955 of file quadtet.cpp.
References mechtop::elements, element::ipp, mechmat::matstiff(), Mm, Mt, mxv(), nne, nodip_quadtet(), and tncomp.
void nod_stresses_ip | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
copies stresses from integration point to nodes
Function computes stresses at nodes of element.
lcid | - load case id | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 894 of file quadtet.cpp.
References vector::a, allocv(), cmulv(), destrv(), mechtop::elements, gauss_points_tet(), mechtop::give_elemnodes(), mechmat::givestress(), intordsm, element::ipp, Mm, Mp, Mt, nne, nod, mechtop::nodes, nodip_quadtet(), node::storestress(), probdesc::stressaver, tncomp, and volumeip().
Referenced by compute_nodestresses(), and stresses().
void node_forces_surf | ( | long | lcid, | |
long | eid, | |||
long * | is, | |||
double * | nv, | |||
vector & | nf | |||
) |
computes nodal forces caused by surface load
Function computes nodal forces caused by surface load.
lcid | - load case id | |
eid | - element id | |
is | - array with surface id | |
nv | - nodal values of surface load | |
nf | - nodal forces due to surface load (output) |
Created by Tomas Koudelka, 11.2008
Definition at line 1804 of file quadtet.cpp.
References matrix::a, ivector::a, vector::a, addv(), allocm(), bf_matrix(), fillm(), fillv(), gauss_points_tr(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), intordb, jac2d_3d(), lgvectortransfblock(), locglob_nodeval(), mechtop::locsystems(), matrix::m, Mt, mxv(), matrix::n, napfun, ndofe, nne, nnj(), nnsurf, nodes, quadtetrahedral_surfnod(), and transf_matrix().
Referenced by loadel::surfaceload().
void nonloc_internal_forces | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci, | |||
vector & | ifor, | |||
vector & | x, | |||
vector & | y, | |||
vector & | z | |||
) |
computes internal forces of the element nonlocally
Function computes nonlocal internal forces.
lcid | - number of load case | |
eid | - element id | |
ri,ci | - row and column indices | |
ifor | - vector of internal forces (output) | |
x,y,z | - vectors of nodal coordinates |
Created by Tomas Koudelka, 11.2008
Definition at line 1150 of file quadtet.cpp.
References compute_nonloc_nlstress(), elem_integration(), and nonlocstress.
Referenced by res_nonloc_internal_forces().
void res_eigstrain_forces | ( | long | lcid, | |
long | eid, | |||
vector & | nfor | |||
) |
computes resulting nodal forces of the element caused by eigenstrains
Function computes resulting contributions from eigenstrains to the right hand side
eid | - element id | |
ifor | - vector of internal forces |
Created by Tomas Koudelka, 11.2008
Definition at line 1338 of file quadtet.cpp.
References copyv(), eigstrain_forces(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), glvectortransf(), mechtop::locsystems(), Mt, ndofe, nne, nodes, and transf_matrix().
Referenced by elem_eigstrain_forces().
void res_incr_internal_forces | ( | long | lcid, | |
long | eid, | |||
vector & | ifor | |||
) |
computes resulting increments of internal forces of the element
Function computes increments of internal forces.
lcid | - load case id | |
eid | - element id | |
ifor | - vector of increments of internal forces (output) |
Created by Tomas Koudelka, 11.2008
Definition at line 1303 of file quadtet.cpp.
References copyv(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), glvectortransf(), incr_internal_forces(), mechtop::locsystems(), Mt, ndofe, nne, nodes, and transf_matrix().
Referenced by elem_incr_internal_forces().
void res_internal_forces | ( | long | lcid, | |
long | eid, | |||
vector & | ifor | |||
) |
computes resulting internal forces of the element
Function computes resulting internal forces (from correct stresses)
lcid | - number of load case | |
eid | - element id | |
ifor | - vector of internal forces |
Created by Tomas Koudelka, 11.2008
Definition at line 1231 of file quadtet.cpp.
References copyv(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), glvectortransf(), internal_forces(), mechtop::locsystems(), Mt, ndofe, nne, nodes, and transf_matrix().
Referenced by elem_internal_forces().
void res_ip_strains | ( | long | lcid, | |
long | eid | |||
) |
computes resulting strains in the integration points
Function computes strains at integration points of given element.
lcid | - load case id | |
eid | - element id |
Created by Tomas Koudelka, 11.2008
Definition at line 512 of file quadtet.cpp.
References allocm(), allocv(), copyv(), destrm(), destrv(), eldispl(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), ip_strains(), lgvectortransf(), mechtop::locsystems(), Mt, ndofe, nne, nodes, and transf_matrix().
Referenced by compute_ipstrains(), and strains().
void res_ip_stresses | ( | long | lcid, | |
long | eid | |||
) |
computes resulting stresses in the integration points
Function computes stresses at integration points of element.
lcid | - load case id | |
eid | - element id |
Created by Tomas Koudelka, 11.2008
Definition at line 810 of file quadtet.cpp.
References ip_stresses().
Referenced by compute_ipstresses(), and stresses().
void res_load_matrix | ( | long | eid, | |
matrix & | lm | |||
) |
computes resulting load matrix
Function computes load matrix. If it is required, nodal values are transformed to the local coordinate systems.
eid | - number of element | |
lm | - load matrix |
Created by Tomas Koudelka, 11.2008
Definition at line 460 of file quadtet.cpp.
References mechtop::give_elemnodes(), glmatrixtransf(), load_matrix(), mechtop::locsystems(), Mt, ndofe, nne, nodes, and transf_matrix().
void res_mass_matrix | ( | long | eid, | |
matrix & | mm | |||
) |
computes resulting mass matrix
Function computes mass matrix without possible transformations at nodes.
eid | - number of element | |
mm | - mass matrix |
Created by Tomas Koudelka, 11.2008
Definition at line 390 of file quadtet.cpp.
References mechtop::give_elemnodes(), glmatrixtransf(), mechtop::locsystems(), mass_matrix(), Mt, ndofe, nne, nodes, and transf_matrix().
void res_nonloc_internal_forces | ( | long | lcid, | |
long | eid, | |||
vector & | ifor | |||
) |
computes resulting internal forces of the element nonlocally
Function computes resulting internal forces for nonlocal models
lcid | - load case id | |
eid | - element id | |
ifor | - vector of internal forces |
Created by Tomas Koudelka, 11.2008
Definition at line 1267 of file quadtet.cpp.
References copyv(), mechtop::give_elemnodes(), mechtop::give_node_coord3d(), glvectortransf(), mechtop::locsystems(), Mt, ndofe, nne, nodes, nonloc_internal_forces(), and transf_matrix().
Referenced by elem_nonloc_internal_forces().
void res_stiffness_matrix | ( | long | eid, | |
matrix & | sm | |||
) |
computes resulting stiffness matrix
Function computes stiffness matrix of one element. If it is required, nodal values are transformed to the local coordinate systems.
eid | - number of element | |
sm | - stiffness matrix (output) |
Created by Tomas Koudelka, 11.2008
Definition at line 316 of file quadtet.cpp.
References mechtop::give_elemnodes(), glmatrixtransf(), mechtop::locsystems(), Mt, ndofe, nne, nodes, stiffness_matrix(), and transf_matrix().
Referenced by stiffmat().
void stiffness_matrix | ( | long | eid, | |
long | ri, | |||
long | ci, | |||
matrix & | sm | |||
) |
computes stiffness matrix
Function computes stiffness matrix of one element.
eid | - number of element | |
ri,ci | - row and column indices | |
sm | - stiffness matrix (output) |
Created by Tomas Koudelka, 11.2008
Definition at line 263 of file quadtet.cpp.
References vector::a, allocv(), bdbjac(), destrv(), mechtop::elements, fillm(), gauss_points_tet(), geom_matrix(), mechtop::give_node_coord3d(), intordsm, element::ipp, mechmat::matstiff(), Mm, Mt, ndofe, nne, and tncomp.
Referenced by res_stiffness_matrix().
void strains | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes required strains
Function computes strains at required positions.
lcid | - load case id | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 732 of file quadtet.cpp.
References allocv(), cenodes, destrv(), enodes, aepoints::give_aepcoord(), aepoints::give_naep(), aepoints::give_ncomp(), aepoints::give_sid(), intpts, Mm, Mp, nne, nod_strains_comp(), nod_strains_ip(), nowhere, print_err(), res_ip_strains(), aepoints::storevalues(), mechmat::stra, probdesc::strainaver, aepoints::tape, tncomp, and userdefined.
void stresses | ( | long | lcid, | |
long | eid, | |||
long | ri, | |||
long | ci | |||
) |
computes required stresses
Function computes stresses at required positions.
lcid | - load case id | |
eid | - element id | |
ri,ci | - row and column indices |
Created by Tomas Koudelka, 11.2008
Definition at line 996 of file quadtet.cpp.
References allocv(), destrv(), enodes, aepoints::give_aepcoord(), aepoints::give_naep(), aepoints::give_ncomp(), aepoints::give_sid(), intpts, Mm, Mp, nod_stresses_ip(), nowhere, print_err(), res_ip_stresses(), aepoints::storevalues(), mechmat::stre, probdesc::stressaver, aepoints::tape, and userdefined.
Referenced by computestresses().
computes transformation matrix
Function assembles transformation matrix from local nodal coordinate system to the global coordinate system x_g = T x_l
nodes | - nodes of element | |
tmat | - transformation matrix (output) |
Created by Tomas Koudelka, 11.2008
Definition at line 220 of file quadtet.cpp.
References node::e1, node::e2, node::e3, fillm(), matrix::m, Mt, ivector::n, mechtop::nodes, and node::transf.
Referenced by nod_strains_comp(), node_forces_surf(), res_eigstrain_forces(), res_incr_internal_forces(), res_internal_forces(), res_ip_strains(), res_load_matrix(), res_mass_matrix(), res_nonloc_internal_forces(), and res_stiffness_matrix().
double volumeip | ( | long | eid, | |
double | w | |||
) |
computes volume appropriate to integration point
Function computes volume appropriate to integration point.
eid | - element id | |
w | - weight of given integration point |
Created by Tomas Koudelka, 3.2009
Definition at line 491 of file quadtet.cpp.
References mechtop::give_volume(), and Mt.
Referenced by nod_eqother_ip(), nod_strains_ip(), and nod_stresses_ip().
long* cncomp |
array containing cumulative numbers of components of stress and strain tensors
Definition at line 137 of file quadtet.h.
Referenced by elem_integration(), ip_strains(), quadtet(), and ~quadtet().
long intordb |
array of numbers of integration points in sets
Definition at line 153 of file quadtet.h.
Referenced by node_forces_surf(), and quadtet().
long intordmm |
order of integration of mass matrix
Definition at line 151 of file quadtet.h.
Referenced by load_matrix(), mass_matrix(), and quadtet().
long** intordsm |
array of orders of integration of stiffness matrix
Definition at line 149 of file quadtet.h.
Referenced by aver_strains(), compute_eigstress(), compute_nlstress(), compute_nlstressincr(), compute_nonloc_nlstress(), elem_integration(), export_gid_gauss_pt(), mechtop::give_intordsm(), inicipval(), intpointval(), intpointval2(), ip_elast_stresses(), ip_strains(), ip_stresses(), ipcoord(), ipvolume(), local_values(), nod_eqother_ip(), nod_strains_ip(), nod_stresses_ip(), quadtet(), stiffness_matrix(), and ~quadtet().
long napfun |
number of approximated functions on the element
Definition at line 139 of file quadtet.h.
Referenced by mechtop::give_napfun(), load_matrix(), mass_matrix(), node_forces_surf(), and quadtet().
long nb |
number of blocks
Definition at line 157 of file quadtet.h.
Referenced by aver_strains(), export_gid_gauss_pt(), mechtop::give_nb(), mechtop::give_nb_te(), inicipval(), intpointval(), intpointval2(), ip_strains(), ipvolume(), quadtet(), and ~quadtet().
long* ncomp |
array containing numbers of components of stress and strain tensors
Definition at line 135 of file quadtet.h.
Referenced by aver_strains(), ip_strains(), quadtet(), and ~quadtet().
long ndofe |
number of DOFs on the element
Definition at line 127 of file quadtet.h.
Referenced by elem_integration(), mechtop::give_ndofe(), ip_elast_stresses(), ip_strains(), load_matrix(), mass_matrix(), nod_strains_comp(), node_forces_surf(), quadtet(), res_eigstrain_forces(), res_incr_internal_forces(), res_internal_forces(), res_ip_strains(), res_load_matrix(), res_mass_matrix(), res_nonloc_internal_forces(), res_stiffness_matrix(), and stiffness_matrix().
long ned |
number of edges
Definition at line 141 of file quadtet.h.
Referenced by mechtop::give_ned(), and quadtet().
long** nip |
array of numbers of integration points on surface
Definition at line 155 of file quadtet.h.
Referenced by mechtop::give_nip(), quadtet(), and ~quadtet().
long nne |
number of nodes on one element
Definition at line 129 of file quadtet.h.
Referenced by approx(), aver_strains(), bf_matrix(), geom_matrix(), mechtop::give_nne(), inicipval(), ip_strains(), ipcoord(), ipvolume(), load_matrix(), mass_matrix(), nod_eqother_ip(), nod_strains_comp(), nod_strains_ip(), nod_stresses_comp(), nod_stresses_ip(), node_forces_surf(), quadtet(), res_eigstrain_forces(), res_incr_internal_forces(), res_internal_forces(), res_ip_strains(), res_load_matrix(), res_mass_matrix(), res_nonloc_internal_forces(), res_stiffness_matrix(), stiffness_matrix(), and strains().
long nned |
number of nodes on one edge
Definition at line 143 of file quadtet.h.
Referenced by mechtop::give_nned(), and quadtet().
long nnsurf |
number of nodes on one surface
Definition at line 147 of file quadtet.h.
Referenced by mechtop::give_nnsurf(), node_forces_surf(), and quadtet().
long nsurf |
number of surfaces
Definition at line 145 of file quadtet.h.
Referenced by mechtop::give_nsurf(), and quadtet().
stress/strain state
Definition at line 159 of file quadtet.h.
Referenced by mechtop::give_ssst(), and quadtet().
long tncomp |
number of components of the strain and stress tensors
Definition at line 131 of file quadtet.h.
Referenced by compute_eigstress(), elem_integration(), mechtop::give_ncomp(), mechtop::give_tncomp(), ip_elast_stresses(), nod_strains_comp(), nod_strains_ip(), nod_stresses_comp(), nod_stresses_ip(), quadtet(), stiffness_matrix(), and strains().
long tnip |
total number of integration points on element
Definition at line 133 of file quadtet.h.
Referenced by mechtop::give_tnip(), and quadtet().