#include <errno.h>
#include "globmat.h"
#include "genfile.h"
#include "mathem.h"
#include "global.h"
#include "elemhead.h"
#include "loadcase.h"
#include "dloadcase.h"
#include "dloadpd.h"
#include "node.h"
#include "element.h"
#include "intpoints.h"
#include "problem.h"
#include "mechmat.h"
#include "creepb.h"
#include "creepbs.h"
#include "mechcrsec.h"
#include "elemswitch.h"
#include "vecttens.h"
Go to the source code of this file.
Functions | |
void | compute_req_val (long lcid) |
void | constr_matrix (long nid, long cid, matrix &lcm) |
void | damping_matrix () |
void | eldispl (long lcid, long eid, double *r) |
void | elprdispl (long lcid, long eid, double *r) |
void | incr_internal_forces (long lcid, double *intfor) |
void | initial_stiffness_matrix (long lcid) |
void | internal_forces (long lcid, double *intfor) |
void | loc_internal_forces (long lcid, double *intfor) |
void | local_global_displ_transf (long lcid) |
void | macrodispl (long lcid, long ncomp, double *r) |
void | mass_matrix (long lcid) |
void | mefel_right_hand_side (long lcid, double *rhs, double *flv) |
void | nodal_eigstrain_forces (long lcid, double *nodfor, double time) |
void | nodal_pore_press_forces (long lcid, double *nodfor, double time) |
void | noddispl (long lcid, double *r, long nid) |
void | nodforce (double *f, long nid, vector &nf) |
void | nonloc_internal_forces (long lcid, double *intfor) |
void | select_noddispl (long lcid, double *r, long nid) |
void | select_nodforce (double *f, long nid, vector &nf) |
void | stiffness_matrix (long lcid) |
void | stress_initdispl (long lcid) |
void compute_req_val | ( | long | lcid | ) |
The function computes all required values this function is called e.g. before output print.
lcid | - load case id |
Created by JK, 1.4.2004 Modified by TKo, 07.2008
Definition at line 1714 of file globmat.cpp.
References compute_ipstrains(), compute_ipstresses(), compute_nodeothers(), compute_nodestrains(), compute_nodestresses(), dloadcase::compute_reactions(), loadcase::compute_reactions(), mechbclc::dlc, earth_pressure, eigen_dynamics, forced_dynamics, growing_mech_structure, mechbclc::lc, lin_floating_subdomain, linear_statics, load_balancing, mat_nonlinear_statics, Mb, mech_timedependent_prob, Mm, Mp, nonlin_floating_subdomain, probdesc::othercomp, probdesc::otherpos, probdesc::otherstate, print_err(), probdesc::reactcomp, probdesc::straincomp, probdesc::strainpos, probdesc::strainstate, probdesc::strcomp, probdesc::stresscomp, probdesc::stresspos, probdesc::stressstate, mechmat::totstrains(), and probdesc::tprob.
Referenced by arclength(), arclengthrv(), arclengthrv1(), difference_method(), displ_control(), garclength(), gnewton_raphson(), gnewton_raphson2(), newmark_method(), newton_raphson(), newton_raphson_coupl_new(), newton_raphson_gparcoupl_lin(), newton_raphson_gparcoupl_nonlin(), newton_raphson_parcoupl_common_dt(), newton_raphson_parcoupl_comp(), newton_raphson_parcoupl_lin(), newton_raphson_parcoupl_nonlin(), newton_raphson_parcoupl_nonlin_new(), newton_raphson_parcoupl_nonlin_old(), newton_raphson_restart(), newton_raphsonrv1(), nonlin_newmark_method(), one_step_arcl(), par_newton_raphson_gparcoupl_lin(), par_newton_raphson_gparcoupl_nonlin(), par_newton_raphson_parcoupl_comp(), par_newton_raphson_parcoupl_lin(), par_newton_raphson_parcoupl_lin_vform(), par_newton_raphson_parcoupl_nonlin(), par_one_step(), par_one_step_mefel(), par_solve_layered_linear_statics(), par_solve_prob_constr_phases(), par_solve_timemech_prob(), par_solve_timemech_prob2(), par_solve_timemech_prob_old(), parallel_solution_eigen_dynamics(), solve_eigen_dynamics(), solve_hemivariational_inequalities(), solve_incremental_floating_subdomains(), solve_linear_floating_subdomains(), solve_linear_floating_subdomains2(), solve_linear_floating_subdomains_old(), solve_linear_statics(), solve_load_balancing(), solve_molecular_dynamics(), solve_prob_constr_phases(), solve_var_stiff_method(), visco_solver(), and visco_solver2().
void constr_matrix | ( | long | nid, | |
long | cid, | |||
matrix & | lcm | |||
) |
The function extracts displacements on one node for given problem.
lcid | - number of load case | |
idn | - node number | |
p | - pointer on problem | |
r | - answer - allocated array for displacement |
Created JK, 11.10.2003 The function assembles constraint matrix of one layer on one node.
nid | - node id | |
cid | - contact id | |
lcm | - local constraint matrix (output) |
Created by JK, 7.12.2002
Definition at line 1492 of file globmat.cpp.
References node::crst, fillm(), mechcrsec::give_onethickness(), Gtm, node::idcs, gtopology::lgnodes, Mc, Mp, Mt, mechtop::nodes, lgnode::nodes, and probdesc::tlm.
Referenced by stiffness_matrix().
void damping_matrix | ( | ) |
The function assembles damping matrix.
Created by JK, 10.8.2005
Definition at line 360 of file globmat.cpp.
References gmatrix::addgm(), probdesc::damp, probdesc::dmass, Dmat, probdesc::dstiff, Gtm, gmatrix::initiate(), massdamp, Mespr, Mmat, Mp, Ndofm, nodamp, Out, gmatrix::prepmat(), print_err(), rayleighdamp, gmatrix::setval(), Smat, probdesc::ssle, stiffdamp, and probdesc::tstorsm.
Referenced by difference_method(), newmark_method(), and nonlin_newmark_method().
void eldispl | ( | long | lcid, | |
long | eid, | |||
double * | r | |||
) |
The function extracts displacements on one element.
lcid | - number of load case | |
eid | - element id | |
r | - allocated array for displacement, it is output parameter |
Created by JK, 9.7.2001
Definition at line 902 of file globmat.cpp.
References vector::a, matrix::a, ivector::a, copyv(), destrv(), mechbclc::dlc, mechtop::elements, forced_dynamics, dloadcase::get_pd(), mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), gtopology::give_ndofn(), gtopology::give_nne(), gtopology::give_node_code_numbers(), gtopology::give_nodes(), growing_mech_structure, Gtm, probdesc::homog, probdesc::lambda, mechbclc::lc, lhsrhs::lhs, Lsrs, mat_nonlinear_statics, Mb, mech_timedependent_prob, Mp, Mt, mxv(), argyrisplate::ndofe, Ndofm, argyrisplate::nne, nod, mechtop::nodedispl, loadcase::pd, reallocv(), element::subtrinitdispl(), probdesc::time, element::tmat, and probdesc::tprob.
Referenced by planeelemqq::compute_error(), planeelemlq::elchar(), linhex::gl_internal_forces(), planeelemlq::gnl_internal_forces(), linhexrot::gnl_internal_forces(), linhex::gnl_internal_forces(), planeelemlq::gnl_stiffness_matrix(), linhexrot::gnl_stiffness_matrix(), linhex::gnl_stiffness_matrix(), mechtop::initial_displ(), linhex::internal_forces(), beamgen3d::internal_forces(), beamel3d::internal_forces(), beamel2d::internal_forces(), beam2dspec::internal_forces(), axisymqq::internal_forces(), axisymlq::internal_forces(), beamgen3d::internal_forces1(), lintet::ip_strains(), planeelemsubqt::mainip_strains(), linhex::mainip_strains(), axisymqq::mainip_strains(), axisymqq::nod_strains(), axisymlq::nod_strains(), quadtet::nod_strains_comp(), quadhex::nod_strains_comp(), planeelemrotlq::nod_strains_comp(), planeelemlt::nod_strains_comp(), planeelemlq::nod_strains_comp(), linhexrot::nod_strains_comp(), linhex::nod_strains_comp(), barelq3d::nod_strains_comp(), barelq2d::nod_strains_comp(), barel2d::nod_strains_comp(), axisymqq::nod_strains_comp(), axisymlt::nod_strains_comp(), axisymqq::nod_strains_old(), axisymqq::nod_stresses(), axisymqq::nod_stresses_comp(), beamgen3d::nodal_displ(), beamel3d::nodal_displ(), beamel2d::nodal_displ(), beamgen3d::nodal_forces(), beamel3d::nodal_forces(), beamel2d::nodal_forces(), planeelemqq::ntdbr_vector(), beamel3d::res_internal_forces(), shelltr::res_ip_strains(), shellq::res_ip_strains(), quadtet::res_ip_strains(), quadhex::res_ip_strains(), q4plate::res_ip_strains(), planeelemrotlt::res_ip_strains(), planeelemrotlq::res_ip_strains(), planeelemqq::res_ip_strains(), planeelemlt::res_ip_strains(), planeelemlq::res_ip_strains(), lintetrot::res_ip_strains(), linhexrot::res_ip_strains(), linhex::res_ip_strains(), dstelem::res_ip_strains(), dktelem::res_ip_strains(), cctelem::res_ip_strains(), barel3d::res_ip_strains(), barel2d::res_ip_strains(), ArgyrisTriangle::res_ip_strains(), quadrilatc::res_mainip_strains(), barelc::res_mainip_strains(), soilplatetr::res_mainip_strains(), soilplateq::res_mainip_strains(), plquadcontact::res_mainip_strains(), planeelemrotlt::res_mainip_strains(), planeelemqt::res_mainip_strains(), beam2dspec::res_mainip_strains(), barelq3d::res_mainip_strains(), barelq2d::res_mainip_strains(), axisymlt::res_mainip_strains(), axisymlq::res_mainip_strains(), beam2dspec::res_mainip_stresses(), mechtop::save_elem_inidispl(), springel::strains(), soilbeam::strains(), beam2dspec::strains(), and beam2dspec::stresses().
void elprdispl | ( | long | lcid, | |
long | eid, | |||
double * | r | |||
) |
The function extracts prescribed displacements on one element.
lcid | - number of load case | |
eid | - element id | |
r | - array of prescribed displacements, it is output parameter |
Created by JK, 9.7.2001
Definition at line 1062 of file globmat.cpp.
References vector::a, matrix::a, ivector::a, copyv(), mechbclc::dlc, mechtop::elements, forced_dynamics, dloadcase::get_pd(), mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), growing_mech_structure, Gtm, probdesc::lambda, mechbclc::lc, mat_nonlinear_statics, Mb, mech_timedependent_prob, Mp, Mt, mxv(), argyrisplate::ndofe, nod, loadcase::pd, reallocv(), element::subtrinitdispl(), probdesc::time, element::tmat, and probdesc::tprob.
Referenced by loadcase::assemble(), and dloadcase::assemble().
void incr_internal_forces | ( | long | lcid, | |
double * | intfor | |||
) |
The function computes vector of increments of internal forces.
lcid | - load case id | |
intfor | - vector of increments of internal forces (output) |
Created by JK, TKo, 29.4.2008
Definition at line 617 of file globmat.cpp.
References vector::a, ivector::a, check_math_errel(), compute_ipstrains(), elem_incr_internal_forces(), mechtop::elements, mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), Gtm, gtopology::leso, locglob(), Mt, mtxv(), argyrisplate::ndofe, Ndofm, mechtop::ne, nullv(), reallocv(), and element::tmat.
Referenced by incr_internal_gforces(), newton_raphson_gparcoupl_lin(), newton_raphson_gparcoupl_nonlin(), newton_raphson_parcoupl_lin(), newton_raphson_parcoupl_nonlin(), one_step(), par_newton_raphson_gparcoupl_lin(), par_newton_raphson_gparcoupl_nonlin(), par_newton_raphson_parcoupl_lin(), par_newton_raphson_parcoupl_lin_vform(), par_newton_raphson_parcoupl_nonlin(), par_one_step(), par_one_step_mefel(), par_solve_prob_constr_phases(), par_solve_timemech_prob(), solve_prob_constr_phases(), and visco_solver().
void initial_stiffness_matrix | ( | long | lcid | ) |
Function assembles initial stiffness matrix.
It is also called geometric stiffness matrix (in R.W. Clough and J. Penzien: Dynamics of structures)
Parameters:
lcid | - load case id |
Created by JK,
Definition at line 303 of file globmat.cpp.
References ivector::a, mechtop::elements, mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), Gtm, gmatrix::initiate(), initstiffmat(), Ismat, gtopology::leso, gmatrix::localize(), Mespr, Mp, Mt, mtxm(), mxm(), argyrisplate::ndofe, Ndofm, mechtop::ne, Out, gmatrix::prepmat(), reallocm(), reallocv(), gmatrix::setval(), probdesc::ssle, element::tmat, and probdesc::tstorsm.
Referenced by solve_linear_stability().
void internal_forces | ( | long | lcid, | |
double * | intfor | |||
) |
The function computes vector of internal forces.
lcid | - load case id | |
intfor | - vetcor of internal forces (output) |
Created by JK, 28.7.2001
Definition at line 406 of file globmat.cpp.
References lin_floating_subdomain, loc_internal_forces(), local, probdesc::matmodel, Mp, nonlin_floating_subdomain, nonloc_internal_forces(), nonlocal, probdesc::otherstate, probdesc::strainstate, probdesc::stressstate, and probdesc::tprob.
Referenced by arclength(), arclength15(), arclengthadapt(), arclengthrv(), arclengthrv1(), loadcase::assemble(), displ_control(), garclength(), general_epressure_varsup(), gnewton_raphson(), gnewton_raphson_one_step(), internal_gforces(), newton_raphson(), newton_raphson_coupl_new(), newton_raphson_gparcoupl_lin(), newton_raphson_gparcoupl_nonlin(), newton_raphson_parcoupl_lin(), newton_raphson_parcoupl_nonlin(), newton_raphson_parcoupl_nonlin_new(), newton_raphson_parcoupl_nonlin_old(), newton_raphson_restart(), newton_raphsonep(), newton_raphsonrv1(), nonlin_newmark_method(), one_step_arcl(), par_gnewton_raphson_one_step(), par_gnewton_raphson_one_step_mefel(), par_newton_raphson_gparcoupl_lin(), par_newton_raphson_gparcoupl_nonlin(), par_newton_raphson_parcoupl_lin(), par_newton_raphson_parcoupl_lin_vform(), par_newton_raphson_parcoupl_nonlin(), par_solve_nonlinear_statics(), par_solve_prob_constr_phases(), par_solve_timemech_prob(), par_solve_timemech_prob_old(), par_visco_mefel_init(), par_visco_solver_init(), solve_linear_statics(), solve_prob_constr_phases(), verlet_method(), visco_solver(), and visco_solver_init().
void loc_internal_forces | ( | long | lcid, | |
double * | intfor | |||
) |
The function computes vector of internal forces, local material models are used.
lcid | - load case id | |
intfor | - vector of internal forces (output) |
Created by JK, TKo, 28.7.2001
Definition at line 448 of file globmat.cpp.
References vector::a, ivector::a, check_math_errel(), compute_ipstrains(), elem_integration_quant(), elem_internal_forces(), mechtop::elements, mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), Gtm, probdesc::homog, gtopology::leso, locglob(), locstress, Mp, Mt, mtxv(), argyrisplate::ncomp, argyrisplate::ndofe, Ndofm, mechtop::ne, nullv(), print_err(), reallocv(), probdesc::strcomp, and element::tmat.
Referenced by internal_forces().
void local_global_displ_transf | ( | long | lcid | ) |
The function transforms nodal displacements form the nodal local coordinate system to the global one.
lcid | - load case id |
Created by JK,
Definition at line 1837 of file globmat.cpp.
References vector::a, allocv(), destrv(), node::e1, node::e2, node::e3, g, mechtop::give_ndofn(), mechtop::give_node_code_numbers(), lhsrhs::lhs, Lsrs, Mt, mxv(), mechtop::nn, noddispl(), mechtop::nodes, and node::transf.
void macrodispl | ( | long | lcid, | |
long | ncomp, | |||
double * | r | |||
) |
The function returns components of displacements vector that are associated with the macro-strains or macro-stresses in case of homogeniztaion problem.
lcid | - load case id | |
ncomp | - number of macro-stress/strain components | |
r | - array of macro-stress/strain components (output parameter) |
Created by Tomas Koudelka, 20.1.2015
Definition at line 1318 of file globmat.cpp.
References copyv(), lhsrhs::give_lhs(), Lsrs, and Ndofm.
void mass_matrix | ( | long | lcid | ) |
Function assembles mass matrix.
Parameters:
lcid | - load case id |
Created by JK,
Definition at line 220 of file globmat.cpp.
References ivector::a, mechtop::elements, mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), Gtm, gmatrix::initiate(), gtopology::leso, gmatrix::localize(), massmat(), Mespr, Mmat, Mp, Mt, mtxm(), mxm(), argyrisplate::ndofe, Ndofm, mechtop::ne, Out, gmatrix::prepmat(), reallocm(), reallocv(), gmatrix::setval(), probdesc::ssle, element::tmat, and probdesc::tstormm.
Referenced by difference_method(), newmark_method(), nonlin_newmark_method(), parallel_solution_eigen_dynamics(), and solve_eigen_dynamics().
void mefel_right_hand_side | ( | long | lcid, | |
double * | rhs, | |||
double * | flv | |||
) |
Function assembles right hand side of system of algebraic equations. It is created by prescribed displacements, forces, moments, eigenstrains, etc. and it contains nodal forces. In case that the flv is not NULL, vector of load caused by forces only is stored in flv
Array rhs has to be cleaned in this function, called functions do not clean the array rhs.
lcid | - load case id | |
rhs | - array of right hand side | |
flv | - array of load vector caused by forces only (default value is NULL) |
Created by JK, Modified by Tomas Koudelka,
Definition at line 1575 of file globmat.cpp.
References dloadcase::assemble(), loadcase::assemble(), mechbclc::dlc, earth_pressure, eigen_dynamics, eigstrain, probdesc::eigstrains, mechmat::est, forced_dynamics, geom_nonlinear_statics, growing_mech_structure, hemivar_inequal, layered_linear_statics, mechbclc::lc, lin_floating_subdomain, linear_stability, linear_statics, load_balancing, mat_nonlinear_statics, Mb, mech_timedependent_prob, Mm, Mp, Ndofm, nodal_eigstrain_forces(), nodal_pore_press_forces(), nonlin_floating_subdomain, nullv(), probdesc::pore_press, print_err(), probdesc::time, and probdesc::tprob.
Referenced by metr_right_hand_side(), newton_raphson_gparcoupl_lin(), newton_raphson_gparcoupl_nonlin(), newton_raphson_parcoupl_lin(), newton_raphson_parcoupl_nonlin(), newton_raphson_parcoupl_nonlin_new(), newton_raphson_parcoupl_nonlin_old(), one_step(), par_newton_raphson_gparcoupl_lin(), par_newton_raphson_gparcoupl_nonlin(), par_newton_raphson_parcoupl_lin(), par_newton_raphson_parcoupl_lin_vform(), par_newton_raphson_parcoupl_nonlin(), par_one_step(), par_one_step_mefel(), par_solve_prob_constr_phases(), par_solve_timemech_prob(), par_solve_timemech_prob_old(), solve_hemivariational_inequalities(), solve_incremental_floating_subdomains(), solve_linear_floating_subdomains(), solve_linear_floating_subdomains2(), solve_linear_floating_subdomains_old(), solve_linear_statics(), solve_load_balancing(), solve_molecular_dynamics(), solve_nonlinear_statics(), solve_prob_constr_phases(), solve_var_stiff_method(), verlet_method(), and visco_solver().
void nodal_eigstrain_forces | ( | long | lcid, | |
double * | nodfor, | |||
double | time | |||
) |
Function adds contributions to internal forces caused by Lagrange multipliers. Function is used in problems with floating subdomains solved by the FETI method.
lcid | - load case id | |
intfor | - array of internal forces (input/output) |
Created by JK, 22.6.2006 The function computes vector of nodal forces caused by eigenstrains. Eigenstrains may contain either eigenstrains or temperature strains.
lcid | - load case id | |
nodfor | - vector of nodal forces (output) | |
time | - actual time |
Created by JK, 28.7.2001
Definition at line 750 of file globmat.cpp.
References vector::a, ivector::a, check_math_errel(), mechbclc::eigstrain_computation(), probdesc::eigstrains, probdesc::eigstrcomp, elem_eigstrain_forces(), mechtop::elements, mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), Gtm, gtopology::leso, locglob(), Mb, Mp, Mt, mtxv(), argyrisplate::ndofe, mechtop::ne, reallocv(), and element::tmat.
Referenced by mefel_right_hand_side(), and loadcase::tempercontrib().
void nodal_pore_press_forces | ( | long | lcid, | |
double * | nodfor, | |||
double | time | |||
) |
The function computes vector of nodal forces caused by eigenstrains. Eigenstrains may contain either eigenstrains or temperature strains.
lcid | - load case id | |
nodfor | - vector of nodal forces (output) | |
time | - actual time |
Created by JK, 28.7.2001
Definition at line 822 of file globmat.cpp.
References vector::a, ivector::a, check_math_errel(), probdesc::eigstrcomp, elem_eigstrain_forces(), mechtop::elements, mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), mechmat::givenonmechq(), Gtm, mechmat::ip, gtopology::leso, locglob(), Mm, Mp, Mt, mtxv(), intpoints::ncompstr, argyrisplate::ndofe, mechtop::ne, pore_pressure, reallocv(), intpoints::ssst, mechmat::storeeigstress(), stress, tensor_vector(), element::tmat, and mechmat::tnip.
Referenced by mefel_right_hand_side().
void noddispl | ( | long | lcid, | |
double * | r, | |||
long | nid | |||
) |
The function extracts displacements on one node.
lcid | - number of load case | |
r | - allocated array for displacement, it is output parameter | |
nid | - node number (node id) |
Created by JK 9.7.2001 Modified by Tomas Koudelka
Definition at line 1181 of file globmat.cpp.
References allocv(), gtopology::approx_weights(), gtopology::give_ndofn(), gtopology::gnodes, Gtm, gnode::masentity, gnode::mnodes, Out, and select_noddispl().
Referenced by general_epressure_varsup(), homogenization(), local_global_displ_transf(), planeelemqq::midpoints(), endnodem::nodal_displacements(), edgem::nodal_displacements(), print_default_dx(), nodeoutm::print_disp(), outdiagm::print_displacements(), mechtop::save_node_inidispl(), mechtop::save_nodval(), write_deflection(), write_displ(), write_gid_displ(), and outdriverm::write_vtk_unkn().
void nodforce | ( | double * | f, | |
long | nid, | |||
vector & | nf | |||
) |
The function extracts nodal force components on one node.
f | - array of components for all nodal forces, it is input parameter | |
nid | - node number (node id) | |
nf | - vector of nodal force, it is output parameter |
Created by Tomas Koudelka 1.10.2013
Definition at line 1339 of file globmat.cpp.
References allocv(), gtopology::approx_weights(), gtopology::give_ndofn(), gtopology::gnodes, Gtm, gnode::masentity, gnode::mnodes, Out, reallocv(), and select_nodforce().
Referenced by write_gid_nforces().
void nonloc_internal_forces | ( | long | lcid, | |
double * | intfor | |||
) |
The function computes vector of internal forces, nonlocal material models are used.
lcid | - load case id | |
intfor | - vector of internal forces (output) |
Created by JK,TKo, 28.7.2001
Definition at line 546 of file globmat.cpp.
References vector::a, ivector::a, check_math_errel(), compute_ipstrains(), elem_local_values(), elem_nonloc_internal_forces(), mechtop::elements, mechtop::give_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), Gtm, gtopology::leso, locglob(), Mm, Mp, Mt, mtxv(), argyrisplate::ndofe, Ndofm, mechtop::ne, mechmat::nonlocaverage(), probdesc::nonlocphase, nullv(), reallocv(), probdesc::strcomp, element::tmat, and mechmat::tnip.
Referenced by internal_forces().
void select_noddispl | ( | long | lcid, | |
double * | r, | |||
long | nid | |||
) |
function selects nodal displacement
lcid | - load case id | |
r | - array of nodal displacements | |
nid | - node id |
20. 8. 2012, JK
Definition at line 1246 of file globmat.cpp.
References mechbclc::dlc, forced_dynamics, dloadcase::get_pd(), mechtop::give_dof(), mechtop::give_ndofn(), growing_mech_structure, probdesc::lambda, mechbclc::lc, lhsrhs::lhs, Lsrs, mat_nonlinear_statics, Mb, mech_timedependent_prob, Mp, Mt, Ndofm, mechtop::nodedispl, loadcase::pd, probdesc::time, and probdesc::tprob.
Referenced by noddispl().
void select_nodforce | ( | double * | f, | |
long | nid, | |||
vector & | nf | |||
) |
The function selects nodal force from the array f and stores its components to the vector nf.
f | - array of components for all nodal forces, it is input parameter | |
nid | - node id | |
nf | - vector of nodal force |
1. 10. 2013, TKo
Definition at line 1399 of file globmat.cpp.
References mechtop::give_dof(), mechtop::give_ndofn(), Mp, Mt, mechtop::nodes, node::r, and probdesc::reactcomp.
Referenced by nodforce(), and mechtop::save_nodforce().
void stiffness_matrix | ( | long | lcid | ) |
Function assembles stiffness matrix.
Parameters:
lcid | - load case id |
Created by JK,
Definition at line 32 of file globmat.cpp.
References ivector::a, lintet::bd_matrix(), check_math_errel(), constr_matrix(), copym(), lintet::dd_matrix(), mechtop::elements, mechtop::give_code_numbers(), mechtop::give_mult_code_numbers(), gtopology::give_ndofe(), mechtop::give_ndofe(), mechtop::give_node_code_numbers(), gmatrix::glocalize(), Gtm, probdesc::homog, gmatrix::initiate(), layered_linear_statics, probdesc::lbtype, gtopology::leso, gtopology::lgnodes, gmatrix::localize(), Ltet, Mespr, Mp, Mt, mtxm(), mxm(), argyrisplate::ndofe, Ndofm, lgnode::ndofn, mechtop::ne, lgnode::nl, mechtop::nln, lgnode::nmult, lgnode::nodes, Out, gmatrix::prepmat(), reallocm(), reallocv(), gmatrix::setval(), Smat, probdesc::ssle, stiffmat(), time, element::tmat, probdesc::tprob, tranm(), and probdesc::tstorsm.
Referenced by arclength(), arclength15(), arclengthadapt(), arclengthrv(), arclengthrv1(), assemble_stiffness_matrix(), difference_method(), displ_control(), newmark_method(), newton_raphson(), newton_raphson_gparcoupl_lin(), newton_raphson_gparcoupl_nonlin(), newton_raphson_parcoupl_lin(), newton_raphson_parcoupl_nonlin(), newton_raphson_parcoupl_nonlin_new(), newton_raphson_parcoupl_nonlin_old(), newton_raphson_restart(), newton_raphsonep(), newton_raphsonrv1(), nonlin_newmark_method(), par_newton_raphson_gparcoupl_lin(), par_newton_raphson_gparcoupl_nonlin(), par_newton_raphson_parcoupl_lin(), par_newton_raphson_parcoupl_lin_vform(), par_newton_raphson_parcoupl_nonlin(), par_solve_layered_linear_statics(), par_solve_linear_statics(), par_solve_nonlinear_statics(), par_solve_prob_constr_phases(), par_solve_timemech_prob(), par_solve_timemech_prob_old(), parallel_solution_eigen_dynamics(), solve_eigen_dynamics(), solve_hemivariational_inequalities(), solve_incremental_floating_subdomains(), solve_layered_linear_statics(), solve_linear_floating_subdomains(), solve_linear_floating_subdomains2(), solve_linear_floating_subdomains_old(), solve_linear_stability(), solve_linear_statics(), solve_load_balancing(), solve_molecular_dynamics(), solve_prob_constr_phases(), solve_var_stiff_method(), and visco_solver().
void stress_initdispl | ( | long | lcid | ) |
The function compute stresses due to initial displacements on the interfaces between new and old part o the structure. It is used due to correct graphics representation of the deformed shape in case of the growing structures.
Parameters:
lcid | - load case id |
Created by TKo, 7.6.2013
Definition at line 1899 of file globmat.cpp.
References compute_ipstrains(), mechtop::elements, mechtop::give_tnip(), mechmat::givestrain(), Gtm, mechmat::ip, element::ipp, gtopology::leso, mechmat::matstiff(), Mm, Mp, Mt, mxv(), argyrisplate::ncomp, intpoints::ncompstr, mechtop::ne, probdesc::otherstate, reallocm(), reallocv(), mechmat::storestress(), probdesc::strainstate, probdesc::stressstate, and argyrisplate::tnip.
Referenced by newton_raphson_gparcoupl_nonlin(), par_newton_raphson_gparcoupl_lin(), par_newton_raphson_gparcoupl_nonlin(), and solve_prob_constr_phases().