📄 elmt_frame2d.c
字号:
/*
* =============================================================================
* ALADDIN Version 1.0 :
* elmt_frame2d.c : Two-dimensional Beam-Column Element
*
* Copyright (C) 1995 by Mark Austin, Xiaoguang Chen, and Wane-Jang Lin
* Institute for Systems Research,
* University of Maryland, College Park, MD 20742
*
* This software is provided "as is" without express or implied warranty.
* Permission is granted to use this software for any on any computer system
* and to redistribute it freely, subject to the following restrictions:
*
* 1. The authors are not responsible for the consequences of use of
* this software, even if they arise from defects in the software.
* 2. The origin of this software must not be misrepresented, either
* by explicit claim or by omission.
* 3. Altered versions must be plainly marked as such, and must not
* be misrepresented as being the original software.
* 4. This notice is to remain intact.
*
* -------------------------------------------------------------------
* Convention for Nodal Forces
* +ve M - anticlockwise(RHT Rule)
* +ve X,Y,Z - along +ve axis
* Convention for Member End Forces
* +ve M - Sagging Moments
* +ve SF - LHS upwards
* +ve AF - Tension(LHS outwards)
* -------------------------------------------------------------------
*
* Written by: Mark Austin, Xiaoguang Chen, and Wane-Jang Lin December 1995
* =============================================================================
*/
#include <math.h>
#include "defs.h"
#include "units.h"
#include "matrix.h"
#include "fe_database.h"
#include "symbol.h"
#include "fe_functions.h"
#include "elmt.h"
/*
#define DEBUG
*/
/* ============================================================== */
/* Element FRAME-2D */
/* 2D Frame Element: Beam_Col Elmt */
/* Frame element : material properties array */
/* Input Properties: */
/* ============================================================== */
/* p->work_material[0] = E;
p->work_material[1] = G;
p->work_material[2] = fy;
p->work_material[3] = ET;
p->work_material[4] = nu;
p->work_material[5] = density;
p->work_material[6] = fu;
p->work_section[0] = Ixx;
p->work_section[1] = Iyy;
p->work_section[2] = Izz;
p->work_section[3] = Ixy;
p->work_section[4] = Ixz;
p->work_section[5] = Iyz;
p->work_section[6] = weight;
p->work_section[7] = bf;
p->work_section[8] = tf;
p->work_section[9] = depth;
p->work_section[10] = area;
p->work_section[11] = plate_thickness;
p->work_section[12] = J;
p->work_section[13] = rT;
p->work_section[14] = width;
p->work_section[15] = tw; */
/* ============================================================== */
#ifdef __STDC__
ARRAY *elmt_frame_2d(ARRAY *p, int isw)
#else
ARRAY *elmt_frame_2d(p, isw)
ARRAY *p;
int isw;
#endif
{
static double nu;
static QUANTITY fy, G, E, ET, density;
static double Ixx, Iyy, Izz, Ixy, Ixz, Iyz, bf, tf, A, depth, weight, EA, EIzz;
double cs, sn, xl,xx,yy, eps,chi,gam, mzc,mass;
double vl1,vl2,tl1,tl2, fx1,fx2,fy1,fy2,mz1,mz2,e6,e2;
double eplas , alp, sig;
int NS_Sub_Incr;
DIMENSIONS *dp_length, *dp_force, *dp_moment;
DIMENSIONS *dp_stress, *dp_degree, *dp_temperature;
double sum, temp;
int i, j, k, ii;
int UNITS_SWITCH;
#ifdef DEBUG
printf("*** Enter elmt_frame_2d() : isw = %4d\n", isw);
#endif
UNITS_SWITCH = CheckUnits();
switch(isw) {
case PROPTY:
/* beam element : material properties */
E.value = p->work_material[0].value;
fy.value = p->work_material[2].value;
ET.value = p->work_material[3].value;
nu = p->work_material[4].value;
density.value = p->work_material[5].value;
if( UNITS_SWITCH == ON ) {
E.dimen = p->work_material[0].dimen;
fy.dimen = p->work_material[2].dimen;
ET.dimen = p->work_material[3].dimen;
density.dimen = p->work_material[5].dimen;
}
/* (1) check poi_ratio value */
if( nu == 0.0 || nu > 0.5 ) {
printf("WARNING >> ... In 2d beam element() - frame_2d - nu value = %9.4f,reset to 0.3 !\n",
nu);
nu = 0.3; /* default poi_ratio value */
}
/* (2) calculate G value */
G.value = p->work_material[1].value = E.value/(1.0 - 2.0*nu) ;
if( UNITS_SWITCH == ON ) G.dimen = E.dimen;
if(E.value/((1.0 - 2.0*nu)) != p->work_material[1].value) {
printf(" elmt_frame_2d(): WARNING: G is not equal to E/(1-2nu), check G for homogeneous material \n");
printf(" : ignore this message for non-homogeneous materials \n");
}
Ixx = p->work_section[0].value;
Iyy = p->work_section[1].value;
Izz = p->work_section[2].value;
Ixy = p->work_section[3].value;
Ixz = p->work_section[4].value;
Iyz = p->work_section[5].value;
weight = p->work_section[6].value;
bf = p->work_section[7].value;
tf = p->work_section[8].value;
depth = p->work_section[9].value;
A = p->work_section[10].value;
EA = E.value*A;
EIzz = E.value*Izz;
break;
case CHERROR:
break;
case STRESS_LOAD:
break;
case STRESS_UPDATE:
break;
case PRESSLD:
break;
case LOAD_MATRIX:
case STRESS:
cs = p->coord[0][1].value - p->coord[0][0].value;
sn = p->coord[1][1].value - p->coord[1][0].value;
xl = sqrt(cs * cs + sn * sn);
cs = cs/xl;
sn = sn/xl;
p->length.value = xl;
xx = 0.5*(p->coord[0][0].value + p->coord[0][1].value); /* xx = 0.5(x1+x2) */
yy = 0.5*(p->coord[1][0].value + p->coord[1][1].value); /* yy = 0.5(y1+y2) */
eps = (cs*(p->displ->uMatrix.daa[0][1] - p->displ->uMatrix.daa[0][0])
+sn*(p->displ->uMatrix.daa[1][1] - p->displ->uMatrix.daa[1][0])) / xl;
/* eps = (u2-u1)cos(phi)+(v2-v1)sin(phi) */
chi = (p->displ->uMatrix.daa[2][1] - p->displ->uMatrix.daa[2][0]) / xl;
gam = -12 *(-sn*(p->displ->uMatrix.daa[0][0] - p->displ->uMatrix.daa[0][1])
+ cs*(p->displ->uMatrix.daa[1][0]- p->displ->uMatrix.daa[1][1]))
/ (xl*xl*xl)-6.*(p->displ->uMatrix.daa[2][0] + p->displ->uMatrix.daa[2][1])/xl/xl;
/* local deflections */
vl1 = -sn * p->displ->uMatrix.daa[0][0] + cs * p->displ->uMatrix.daa[1][0];
vl2 = -sn * p->displ->uMatrix.daa[0][1] + cs * p->displ->uMatrix.daa[1][1];
tl1 = p->displ->uMatrix.daa[2][0];
tl2 = p->displ->uMatrix.daa[2][1];
/* computer axial forces fx, transverse forces fy, and moment mz */
e6 = 6* EIzz /xl/xl;
e2 = 2* EIzz /xl;
/* Elasto_Plastic load case */
/* Considering uniaxial element with x_displ, Fx force */
fx1 = - EA * eps;
fy1 = 2*e6/xl*(vl1-vl2) + e6*(tl1 + tl2);
fx2 = - fx1;
fy2 = - fy1;
mz1 = e6 * ( vl1 - vl2) + e2 * ( 2 * tl1 + tl2 );
mz2 = e6 * ( vl1 - vl2) + e2 * ( tl1 + 2 * tl2 );
mzc = (mz1 - mz2)/2;
/* Add FEF if elmt loaded */
if( p->elmt_load_ptr != NULL) {
p = sld07(p, STRESS);
/* Add FEF to joint forces */
fx1 = fx1 - p->nodal_loads[0].value;
fy1 = fy1 - p->nodal_loads[1].value;
mz1 = mz1 - p->nodal_loads[2].value;
fx2 = fx2 - p->nodal_loads[3].value;
fy2 = fy2 - p->nodal_loads[4].value;
mz2 = mz2 - p->nodal_loads[5].value;
}
/* Assign force values */
if( UNITS_SWITCH == ON ) {
if(UNITS_TYPE == SI) {
dp_length = DefaultUnits("m");
dp_force = DefaultUnits("N");
}
else if(UNITS_TYPE == US) {
dp_length = DefaultUnits("in");
dp_force = DefaultUnits("lbf");
}
/* node no 1 */
UnitsCopy( p->nodal_loads[0].dimen, dp_force );
UnitsCopy( p->nodal_loads[1].dimen, dp_force );
UnitsMultRep( p->nodal_loads[2].dimen, dp_force, dp_length );
/* node no = 2 */
UnitsCopy( p->nodal_loads[3].dimen, p->nodal_loads[0].dimen );
UnitsCopy( p->nodal_loads[4].dimen, p->nodal_loads[1].dimen );
UnitsCopy( p->nodal_loads[5].dimen, p->nodal_loads[2].dimen );
}
p->nodal_loads[0].value = fx1;
p->nodal_loads[1].value = fy1;
p->nodal_loads[2].value = mz1;
p->nodal_loads[3].value = fx2;
p->nodal_loads[4].value = fy2;
p->nodal_loads[5].value = mz2;
if(isw == LOAD_MATRIX ) {
p->nodal_loads[0].value = fx1*cs - fy1*sn;
p->nodal_loads[1].value = fx1*sn + fy1*cs;
p->nodal_loads[2].value = mz1;
p->nodal_loads[3].value = fx2*cs - fy2*sn;
p->nodal_loads[4].value = fx2*sn + fy2*cs;
p->nodal_loads[5].value = mz2;
}
if(isw == STRESS && PRINT_STRESS == ON) {
printf("\n");
printf("Elmt No %3d : \n", p->elmt_no);
if( UNITS_SWITCH == ON ) {
printf("Coords (X,Y) = (%10.3f %s, %10.3f %s)\n",
xx/dp_length->scale_factor,dp_length->units_name,
yy/dp_length->scale_factor,dp_length->units_name);
printf("exx = %13.5e , curva = %13.5g , gamma = %13.5e\n", eps,chi, gam);
printf("\n");
/* node i */
printf(" Fx1 = %13.5e %s Fy1 = %13.5e %s Mz1 = %13.5e %s\n",
p->nodal_loads[0].value/p->nodal_loads[0].dimen->scale_factor,
p->nodal_loads[0].dimen->units_name,
p->nodal_loads[1].value/p->nodal_loads[1].dimen->scale_factor,
p->nodal_loads[1].dimen->units_name,
p->nodal_loads[2].value/p->nodal_loads[2].dimen->scale_factor,
p->nodal_loads[2].dimen->units_name);
/* node j */
printf(" Fx2 = %13.5e %s Fy2 = %13.5e %s Mz2 = %13.5e %s\n",
p->nodal_loads[3].value/p->nodal_loads[3].dimen->scale_factor,
p->nodal_loads[3].dimen->units_name,
p->nodal_loads[4].value/p->nodal_loads[4].dimen->scale_factor,
p->nodal_loads[4].dimen->units_name,
p->nodal_loads[5].value/p->nodal_loads[5].dimen->scale_factor,
p->nodal_loads[5].dimen->units_name);
printf("\n");
/* Member Forces */
printf(" Axial Force : x-direction = %13.5e %s\n",
-p->nodal_loads[0].value/p->nodal_loads[0].dimen->scale_factor,
p->nodal_loads[0].dimen->units_name);
printf(" Shear Force : y-direction = %13.5e %s\n",
p->nodal_loads[1].value/p->nodal_loads[1].dimen->scale_factor,
p->nodal_loads[1].dimen->units_name);
printf("\n");
free((char *) dp_length->units_name);
free((char *) dp_length);
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -